Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Ifpack2_Details_TriDiSolver_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_DETAILS_TRIDISOLVER_DEF_HPP
44 #define IFPACK2_DETAILS_TRIDISOLVER_DEF_HPP
45 
46 #include "Ifpack2_LocalFilter.hpp"
47 #include "Teuchos_LAPACK.hpp"
48 
49 #ifdef HAVE_MPI
50 # include <mpi.h>
52 #else
53 # include "Teuchos_DefaultSerialComm.hpp"
54 #endif // HAVE_MPI
55 
56 
57 namespace Ifpack2 {
58 namespace Details {
59 
61 // Non-stub (full) implementation
63 
64 template<class MatrixType>
67  A_ (A),
68  initializeTime_ (0.0),
69  computeTime_ (0.0),
70  applyTime_ (0.0),
71  numInitialize_ (0),
72  numCompute_ (0),
73  numApply_ (0),
74  isInitialized_ (false),
75  isComputed_ (false)
76 {}
77 
78 
79 template<class MatrixType>
83  A_.is_null (), std::runtime_error, "Ifpack2::Details::TriDiSolver::"
84  "getDomainMap: The input matrix A is null. Please call setMatrix() with a "
85  "nonnull input matrix before calling this method.");
86  // For an input matrix A, TriDiSolver solves Ax=b for x.
87  // Thus, its Maps are reversed from those of the input matrix.
88  return A_->getRangeMap ();
89 }
90 
91 
92 template<class MatrixType>
96  A_.is_null (), std::runtime_error, "Ifpack2::Details::TriDiSolver::"
97  "getRangeMap: The input matrix A is null. Please call setMatrix() with a "
98  "nonnull input matrix before calling this method.");
99  // For an input matrix A, TriDiSolver solves Ax=b for x.
100  // Thus, its Maps are reversed from those of the input matrix.
101  return A_->getDomainMap ();
102 }
103 
104 
105 template<class MatrixType>
106 void
109  (void) params; // this preconditioner doesn't currently take any parameters
110 }
111 
112 
113 template<class MatrixType>
114 bool
116  return isInitialized_;
117 }
118 
119 
120 template<class MatrixType>
121 bool
123  return isComputed_;
124 }
125 
126 
127 template<class MatrixType>
128 int
130  return numInitialize_;
131 }
132 
133 
134 template<class MatrixType>
135 int
137  return numCompute_;
138 }
139 
140 
141 template<class MatrixType>
142 int
144  return numApply_;
145 }
146 
147 
148 template<class MatrixType>
149 double
151  return initializeTime_;
152 }
153 
154 
155 template<class MatrixType>
156 double
158  return computeTime_;
159 }
160 
161 
162 template<class MatrixType>
163 double
165  return applyTime_;
166 }
167 
168 
169 template<class MatrixType>
172  return A_;
173 }
174 
175 
176 template<class MatrixType>
178 reset ()
179 {
180  isInitialized_ = false;
181  isComputed_ = false;
182  A_local_ = Teuchos::null;
183  A_local_tridi_.reshape (0);
184  ipiv_.resize (0);
185 }
186 
187 
188 template<class MatrixType>
191 {
192  // It's legitimate to call setMatrix() with a null input. This has
193  // the effect of resetting the preconditioner's internal state.
194  if (! A_.is_null ()) {
195  const global_size_t numRows = A->getRangeMap ()->getGlobalNumElements ();
196  const global_size_t numCols = A->getDomainMap ()->getGlobalNumElements ();
198  numRows != numCols, std::invalid_argument, "Ifpack2::Details::TriDiSolver::"
199  "setMatrix: Input matrix must be (globally) square. "
200  "The matrix you provided is " << numRows << " by " << numCols << ".");
201  }
202  // Clear any previously computed objects.
203  reset ();
204 
205  // Now that we've cleared the state, we can keep the matrix.
206  A_ = A;
207 }
208 
209 
210 template<class MatrixType>
212 {
213  using Teuchos::Comm;
214  using Teuchos::null;
215  using Teuchos::RCP;
216  using Teuchos::rcp;
217  using Teuchos::Time;
218  using Teuchos::TimeMonitor;
219  const std::string timerName ("Ifpack2::Details::TriDiSolver::initialize");
220 
221  RCP<Time> timer = TimeMonitor::lookupCounter (timerName);
222  if (timer.is_null ()) {
223  timer = TimeMonitor::getNewCounter (timerName);
224  }
225 
226  { // Begin timing here.
227  Teuchos::TimeMonitor timeMon (*timer);
228 
230  A_.is_null (), std::runtime_error, "Ifpack2::Details::TriDiSolver::"
231  "initialize: The input matrix A is null. Please call setMatrix() "
232  "with a nonnull input before calling this method.");
233 
235  ! A_->hasColMap (), std::invalid_argument, "Ifpack2::Details::TriDiSolver: "
236  "The constructor's input matrix must have a column Map, "
237  "so that it has local indices.");
238 
239  // Clear any previously computed objects.
240  reset ();
241 
242  // Make the local filter of the input matrix A.
243  if (A_->getComm ()->getSize () > 1) {
244  A_local_ = rcp (new LocalFilter<row_matrix_type> (A_));
245  } else {
246  A_local_ = A_;
247  }
248 
250  A_local_.is_null (), std::logic_error, "Ifpack2::Details::TriDiSolver::"
251  "initialize: A_local_ is null after it was supposed to have been "
252  "initialized. Please report this bug to the Ifpack2 developers.");
253 
254  // Allocate the TriDi local matrix and the pivot array.
255  const size_t numRows = A_local_->getNodeNumRows ();
256  const size_t numCols = A_local_->getNodeNumCols ();
258  numRows != numCols, std::logic_error, "Ifpack2::Details::TriDiSolver::"
259  "initialize: Local filter matrix is not square. This should never happen. "
260  "Please report this bug to the Ifpack2 developers.");
261  A_local_tridi_.reshape (numRows);
262  ipiv_.resize (numRows);
263  std::fill (ipiv_.begin (), ipiv_.end (), 0);
264 
265  isInitialized_ = true;
266  ++numInitialize_;
267  }
268 
269  // timer->totalElapsedTime() returns the total time over all timer
270  // calls. Thus, we use = instead of +=.
271  initializeTime_ = timer->totalElapsedTime ();
272 }
273 
274 
275 template<class MatrixType>
277 {}
278 
279 
280 template<class MatrixType>
282 {
283  using Teuchos::RCP;
284  const std::string timerName ("Ifpack2::Details::TriDiSolver::compute");
285 
286  RCP<Teuchos::Time> timer = Teuchos::TimeMonitor::lookupCounter (timerName);
287  if (timer.is_null ()) {
288  timer = Teuchos::TimeMonitor::getNewCounter (timerName);
289  }
290 
291  // Begin timing here.
292  {
293  Teuchos::TimeMonitor timeMon (*timer);
295  A_.is_null (), std::runtime_error, "Ifpack2::Details::TriDiSolver::"
296  "compute: The input matrix A is null. Please call setMatrix() with a "
297  "nonnull input, then call initialize(), before calling this method.");
298 
300  A_local_.is_null (), std::logic_error, "Ifpack2::Details::TriDiSolver::"
301  "compute: A_local_ is null. Please report this bug to the Ifpack2 "
302  "developers.");
303 
304  isComputed_ = false;
305  if (! this->isInitialized ()) {
306  this->initialize ();
307  }
308  extract (A_local_tridi_, *A_local_); // extract the tridi local matrix
309 
310  factor (A_local_tridi_, ipiv_ ()); // factor the tridi local matrix
311 
312  isComputed_ = true;
313  ++numCompute_;
314  }
315  // timer->totalElapsedTime() returns the total time over all timer
316  // calls. Thus, we use = instead of +=.
317  computeTime_ = timer->totalElapsedTime ();
318 }
319 
320 template<class MatrixType>
322  const Teuchos::ArrayView<int>& ipiv)
323 {
324  // Fill the LU permutation array with zeros.
325  std::fill (ipiv.begin (), ipiv.end (), 0);
326 
328  int INFO = 0;
329  lapack.GTTRF (A.numRowsCols (),
330  A.DL(),
331  A.D(),
332  A.DU(),
333  A.DU2(),
334  ipiv.getRawPtr (), &INFO);
335  // INFO < 0 is a bug.
337  INFO < 0, std::logic_error, "Ifpack2::Details::TriDiSolver::factor: "
338  "LAPACK's _GTTRF (tridiagonal LU factorization with partial pivoting) "
339  "was called incorrectly. INFO = " << INFO << " < 0. "
340  "Please report this bug to the Ifpack2 developers.");
341  // INFO > 0 means the matrix is singular. This is probably an issue
342  // either with the choice of rows the rows we extracted, or with the
343  // input matrix itself.
345  INFO > 0, std::runtime_error, "Ifpack2::Details::TriDiSolver::factor: "
346  "LAPACK's _GTTRF (tridiagonal LU factorization with partial pivoting) "
347  "reports that the computed U factor is exactly singular. U(" << INFO <<
348  "," << INFO << ") (one-based index i) is exactly zero. This probably "
349  "means that the input matrix has a singular diagonal block.");
350 }
351 
352 
353 template<class MatrixType>
354 void TriDiSolver<MatrixType, false>::
355 applyImpl (const MV& X,
356  MV& Y,
357  const Teuchos::ETransp mode,
358  const scalar_type alpha,
359  const scalar_type beta) const
360 {
361  using Teuchos::ArrayRCP;
362  using Teuchos::RCP;
363  using Teuchos::rcp;
364  using Teuchos::rcpFromRef;
365  using Teuchos::CONJ_TRANS;
366  using Teuchos::TRANS;
367 
368  const int numVecs = static_cast<int> (X.getNumVectors ());
369  if (alpha == STS::zero ()) { // don't need to solve the linear system
370  if (beta == STS::zero ()) {
371  // Use BLAS AXPY semantics for beta == 0: overwrite, clobbering
372  // any Inf or NaN values in Y (rather than multiplying them by
373  // zero, resulting in NaN values).
374  Y.putScalar (STS::zero ());
375  }
376  else { // beta != 0
377  Y.scale (STS::zero ());
378  }
379  }
380  else { // alpha != 0; must solve the linear system
382  // If beta is nonzero, Y is not constant stride, or alpha != 1, we
383  // have to use a temporary output multivector Y_tmp. It gets a
384  // copy of alpha*X, since GETRS overwrites its (multi)vector input
385  // with its output.
386  RCP<MV> Y_tmp;
387  if (beta == STS::zero () && Y.isConstantStride () && alpha == STS::one ()) {
388  deep_copy(Y, X);
389  Y_tmp = rcpFromRef (Y);
390  }
391  else {
392  Y_tmp = rcp (new MV (createCopy(X))); // constructor copies X
393  if (alpha != STS::one ()) {
394  Y_tmp->scale (alpha);
395  }
396  }
397  const int Y_stride = static_cast<int> (Y_tmp->getStride ());
398  ArrayRCP<scalar_type> Y_view = Y_tmp->get1dViewNonConst ();
399  scalar_type* const Y_ptr = Y_view.getRawPtr ();
400  int INFO = 0;
401  const char trans =
402  (mode == CONJ_TRANS ? 'C' : (mode == TRANS ? 'T' : 'N'));
403  lapack.GTTRS (trans, A_local_tridi_.numRowsCols(), numVecs,
404  A_local_tridi_.DL(),
405  A_local_tridi_.D(),
406  A_local_tridi_.DU(),
407  A_local_tridi_.DU2(),
408  ipiv_.getRawPtr (), Y_ptr, Y_stride, &INFO);
410  INFO != 0, std::runtime_error, "Ifpack2::Details::TriDiSolver::"
411  "applyImpl: LAPACK's _GTTRS (tridiagonal solve using LU factorization "
412  "with partial pivoting) failed with INFO = " << INFO << " != 0.");
413 
414  if (beta != STS::zero ()) {
415  Y.update (alpha, *Y_tmp, beta);
416  }
417  else if (! Y.isConstantStride ()) {
418  deep_copy(Y, *Y_tmp);
419  }
420  }
421 }
422 
423 
424 template<class MatrixType>
426 apply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
427  Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
428  Teuchos::ETransp mode,
429  scalar_type alpha,
430  scalar_type beta) const
431 {
432  using Teuchos::ArrayView;
433  using Teuchos::as;
434  using Teuchos::RCP;
435  using Teuchos::rcp;
436  using Teuchos::rcpFromRef;
437 
438  const std::string timerName ("Ifpack2::Details::TriDiSolver::apply");
439  RCP<Teuchos::Time> timer = Teuchos::TimeMonitor::lookupCounter (timerName);
440  if (timer.is_null ()) {
441  timer = Teuchos::TimeMonitor::getNewCounter (timerName);
442  }
443 
444  // Begin timing here.
445  {
446  Teuchos::TimeMonitor timeMon (*timer);
447 
449  ! isComputed_, std::runtime_error, "Ifpack2::Details::TriDiSolver::apply: "
450  "You must have called the compute() method before you may call apply(). "
451  "You may call the apply() method as many times as you want after calling "
452  "compute() once, but you must have called compute() at least once.");
453 
454  const size_t numVecs = X.getNumVectors ();
455 
457  numVecs != Y.getNumVectors (), std::runtime_error,
458  "Ifpack2::Details::TriDiSolver::apply: X and Y have different numbers "
459  "of vectors. X has " << X.getNumVectors () << ", but Y has "
460  << X.getNumVectors () << ".");
461 
462  if (numVecs == 0) {
463  return; // done! nothing to do
464  }
465 
466  // Set up "local" views of X and Y.
467  RCP<const MV> X_local;
468  RCP<MV> Y_local;
469  const bool multipleProcs = (A_->getRowMap ()->getComm ()->getSize () >= 1);
470  if (multipleProcs) {
471  // Interpret X and Y as "local" multivectors, that is, in the
472  // local filter's domain resp. range Maps. "Interpret" means that
473  // we create views with different Maps; we don't have to copy.
474  X_local = X.offsetView (A_local_->getDomainMap (), 0);
475  Y_local = Y.offsetViewNonConst (A_local_->getRangeMap (), 0);
476  }
477  else { // only one process in A_'s communicator
478  // X and Y are already "local"; no need to set up local views.
479  X_local = rcpFromRef (X);
480  Y_local = rcpFromRef (Y);
481  }
482 
483  // Apply the local operator:
484  // Y_local := beta*Y_local + alpha*M^{-1}*X_local
485  this->applyImpl (*X_local, *Y_local, mode, alpha, beta);
486 
487  ++numApply_; // We've successfully finished the work of apply().
488  }
489 
490  // timer->totalElapsedTime() returns the total time over all timer
491  // calls. Thus, we use = instead of +=.
492  applyTime_ = timer->totalElapsedTime ();
493 }
494 
495 
496 template<class MatrixType>
498 {
499  std::ostringstream out;
500 
501  // Output is a valid YAML dictionary in flow style. If you don't
502  // like everything on a single line, you should call describe()
503  // instead.
504  out << "\"Ifpack2::Details::TriDiSolver\": ";
505  out << "{";
506  if (this->getObjectLabel () != "") {
507  out << "Label: \"" << this->getObjectLabel () << "\", ";
508  }
509  out << "Initialized: " << (isInitialized () ? "true" : "false") << ", "
510  << "Computed: " << (isComputed () ? "true" : "false") << ", ";
511 
512  if (A_.is_null ()) {
513  out << "Matrix: null";
514  }
515  else {
516  out << "Matrix: not null"
517  << ", Global matrix dimensions: ["
518  << A_->getGlobalNumRows () << ", " << A_->getGlobalNumCols () << "]";
519  }
520 
521  out << "}";
522  return out.str ();
523 }
524 
525 
526 template<class MatrixType>
528  const Teuchos::EVerbosityLevel verbLevel) const {
529  using Teuchos::FancyOStream;
530  using Teuchos::OSTab;
531  using Teuchos::RCP;
532  using Teuchos::rcpFromRef;
533  using std::endl;
534 
535  if (verbLevel == Teuchos::VERB_NONE) {
536  return;
537  }
538  else {
539  RCP<FancyOStream> ptrOut = rcpFromRef (out);
540  OSTab tab1 (ptrOut);
541  if (this->getObjectLabel () != "") {
542  out << "label: " << this->getObjectLabel () << endl;
543  }
544  out << "initialized: " << (isInitialized_ ? "true" : "false") << endl
545  << "computed: " << (isComputed_ ? "true" : "false") << endl
546  << "number of initialize calls: " << numInitialize_ << endl
547  << "number of compute calls: " << numCompute_ << endl
548  << "number of apply calls: " << numApply_ << endl
549  << "total time in seconds in initialize: " << initializeTime_ << endl
550  << "total time in seconds in compute: " << computeTime_ << endl
551  << "total time in seconds in apply: " << applyTime_ << endl;
552  if (verbLevel >= Teuchos::VERB_EXTREME) {
553  out << "A_local_tridi_:" << endl;
554  A_local_tridi_.print(out);
555  }
556  out << "ipiv_: " << Teuchos::toString (ipiv_) << endl;
557  }
558 }
559 
560 template<class MatrixType>
562  const Teuchos::EVerbosityLevel verbLevel) const
563 {
564  using Teuchos::FancyOStream;
565  using Teuchos::OSTab;
566  using Teuchos::RCP;
567  using Teuchos::rcpFromRef;
568  using std::endl;
569 
570  RCP<FancyOStream> ptrOut = rcpFromRef (out);
571  OSTab tab0 (ptrOut);
572  if (A_.is_null ()) {
573  // If A_ is null, we don't have a communicator, so we can't
574  // safely print local data on all processes. Just print the
575  // local data without arbitration between processes, and hope
576  // for the best.
577  if (verbLevel > Teuchos::VERB_NONE) {
578  out << "Ifpack2::Details::TriDiSolver:" << endl;
579  }
580  describeLocal (out, verbLevel);
581  }
582  else {
583  // If A_ is not null, we have a communicator, so we can
584  // arbitrate among all processes to print local data.
585  const Teuchos::Comm<int>& comm = * (A_->getRowMap ()->getComm ());
586  const int myRank = comm.getRank ();
587  const int numProcs = comm.getSize ();
588  if (verbLevel > Teuchos::VERB_NONE && myRank == 0) {
589  out << "Ifpack2::Details::TriDiSolver:" << endl;
590  }
591  OSTab tab1 (ptrOut);
592  for (int p = 0; p < numProcs; ++p) {
593  if (myRank == p) {
594  out << "Process " << myRank << ":" << endl;
595  describeLocal (out, verbLevel);
596  }
597  comm.barrier ();
598  comm.barrier ();
599  comm.barrier ();
600  } // for p = 0 .. numProcs-1
601  }
602 }
603 
604 template<class MatrixType>
606  const row_matrix_type& A_local)
607 {
608  using Teuchos::Array;
609  using Teuchos::ArrayView;
610  typedef local_ordinal_type LO;
611  typedef typename Teuchos::ArrayView<LO>::size_type size_type;
612 
613  // Fill the local tridi matrix with zeros.
614  A_local_tridi.putScalar (STS::zero ());
615 
616  //
617  // Map both row and column indices to local indices. We can use the
618  // row Map's local indices for row indices, and the column Map's
619  // local indices for column indices. It doesn't really matter;
620  // either way is just a permutation of rows and columns.
621  //
622  const map_type& rowMap = * (A_local.getRowMap ());
623 
624  // Temporary arrays to hold the indices and values of the entries in
625  // each row of A_local.
626  const size_type maxNumRowEntries =
627  static_cast<size_type> (A_local.getNodeMaxNumRowEntries ());
628  Array<LO> localIndices (maxNumRowEntries);
629  Array<scalar_type> values (maxNumRowEntries);
630 
631  const LO numLocalRows = static_cast<LO> (rowMap.getNodeNumElements ());
632  const LO minLocalRow = rowMap.getMinLocalIndex ();
633  // This slight complication of computing the upper loop bound avoids
634  // issues if the row Map has zero entries on the calling process.
635  const LO maxLocalRow = minLocalRow + numLocalRows; // exclusive bound
636  for (LO localRow = minLocalRow; localRow < maxLocalRow; ++localRow) {
637  // The LocalFilter automatically excludes "off-process" entries.
638  // That means all the column indices in this row belong to the
639  // domain Map. We can, therefore, just use the local row and
640  // column indices to put each entry directly in the tridi matrix.
641  // It's OK if the column Map puts the local indices in a different
642  // order; the Import will bring them into the correct order.
643  const size_type numEntriesInRow =
644  static_cast<size_type> (A_local.getNumEntriesInLocalRow (localRow));
645  size_t numEntriesOut = 0; // ignored
646  A_local.getLocalRowCopy (localRow,
647  localIndices (0, numEntriesInRow),
648  values (0, numEntriesInRow),
649  numEntriesOut);
650  for (LO k = 0; k < numEntriesInRow; ++k) {
651  const LO localCol = localIndices[k];
652  const scalar_type val = values[k];
653  // We use += instead of =, in case there are duplicate entries
654  // in the row. There should not be, but why not be general?
655  // NOTE: we only extract the TriDi part of the row matrix. Do not extract DU2
656  if( localCol >= localRow-1 && localCol <= localRow+1 )
657  A_local_tridi(localRow, localCol) += val;
658  }
659  }
660 }
661 
663 // Stub implementation
665 
666 template<class MatrixType>
667 TriDiSolver<MatrixType, true>::TriDiSolver (const Teuchos::RCP<const row_matrix_type>& A) {
668  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
669 }
670 
671 
672 template<class MatrixType>
675  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
676 }
677 
678 
679 template<class MatrixType>
682  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
683 }
684 
685 
686 template<class MatrixType>
687 void
689  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
690 }
691 
692 
693 template<class MatrixType>
694 bool
696  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
697 }
698 
699 
700 template<class MatrixType>
701 bool
703  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
704 }
705 
706 
707 template<class MatrixType>
708 int
710  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
711 }
712 
713 
714 template<class MatrixType>
715 int
717  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
718 }
719 
720 
721 template<class MatrixType>
722 int
724  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
725 }
726 
727 
728 template<class MatrixType>
729 double
731  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
732 }
733 
734 
735 template<class MatrixType>
736 double
738  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
739 }
740 
741 
742 template<class MatrixType>
743 double
745  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
746 }
747 
748 
749 template<class MatrixType>
752  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
753 }
754 
755 
756 template<class MatrixType>
758 {
759  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
760 }
761 
762 
763 template<class MatrixType>
765 {
766  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
767 }
768 
769 
770 template<class MatrixType>
771 TriDiSolver<MatrixType, true>::~TriDiSolver ()
772 {
773  // Destructors should never throw exceptions.
774 }
775 
776 
777 template<class MatrixType>
779 {
780  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
781 }
782 
783 
784 template<class MatrixType>
786  const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
787  Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
788  Teuchos::ETransp mode,
789  scalar_type alpha,
790  scalar_type beta) const
791 {
792  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
793 }
794 
795 
796 template<class MatrixType>
797 std::string
798 TriDiSolver<MatrixType, true>::description () const
799 {
800  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
801 }
802 
803 
804 template<class MatrixType>
805 void TriDiSolver<MatrixType, true>::describe(Teuchos::FancyOStream& out,
806  const Teuchos::EVerbosityLevel verbLevel) const
807 {
808  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Not implemented");
809 }
810 
811 }// namespace Details
812 } // namespace Ifpack2
813 
814 #define IFPACK2_DETAILS_TRIDISOLVER_INSTANT(S,LO,GO,N) \
815  template class Ifpack2::Details::TriDiSolver< Tpetra::RowMatrix<S, LO, GO, N> >;
816 
817 #endif // IFPACK2_DETAILS_TRIDISOLVER_HPP
virtual Teuchos::RCP< const Tpetra::Map< MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > getDomainMap() const =0
The domain Map of this operator.
virtual int getSize() const =0
virtual Teuchos::RCP< const Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > getMatrix() const =0
The input matrix given to the constructor.
OrdinalType numRowsCols() const
virtual int getRank() const =0
&quot;Preconditioner&quot; that uses LAPACK&#39;s tridi LU.
Definition: Ifpack2_Details_TriDiSolver_decl.hpp:74
basic_OSTab< char > OSTab
iterator begin() const
static RCP< Time > getNewCounter(const std::string &name)
static RCP< Time > lookupCounter(const std::string &name)
void GTTRS(const char &TRANS, const OrdinalType &n, const OrdinalType &nrhs, const ScalarType *dl, const ScalarType *d, const ScalarType *du, const ScalarType *du2, const OrdinalType *IPIV, ScalarType *B, const OrdinalType &ldb, OrdinalType *info) const
void GTTRF(const OrdinalType &n, ScalarType *dl, ScalarType *d, ScalarType *du, ScalarType *du2, OrdinalType *IPIV, OrdinalType *info) const
basic_FancyOStream< char > FancyOStream
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
MatrixType::scalar_type scalar_type
The type of entries in the input (global) matrix.
Definition: Ifpack2_Details_TriDiSolver_decl.hpp:107
virtual Teuchos::RCP< const Tpetra::Map< MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > getRangeMap() const =0
The range Map of this operator.
virtual void barrier() const =0
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
TriDiSolver(const Teuchos::RCP< const row_matrix_type > &matrix)
Constructor.
Definition: Ifpack2_Details_TriDiSolver_def.hpp:66
T * getRawPtr() const
virtual bool isComputed() const =0
True if the preconditioner has been successfully computed, else false.
virtual void setParameters(const Teuchos::ParameterList &List)=0
Set this preconditioner&#39;s parameters.
iterator end() const
TypeTo as(const TypeFrom &t)
int putScalar(const ScalarType value=Teuchos::ScalarTraits< ScalarType >::zero())
Access only local rows and columns of a sparse matrix.
Definition: Ifpack2_LocalFilter_decl.hpp:160
virtual void apply(const Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &X, Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, MatrixType::scalar_typealpha=Teuchos::ScalarTraits< MatrixType::scalar_type >::one(), MatrixType::scalar_typebeta=Teuchos::ScalarTraits< MatrixType::scalar_type >::zero()) const =0
Apply the preconditioner to X, putting the result in Y.
virtual bool isInitialized() const =0
True if the preconditioner has been successfully initialized, else false.
std::string toString(const T &t)
virtual void setMatrix(const Teuchos::RCP< const Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > &A)=0
Set the new matrix.