Ifpack Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Ifpack_AdditiveSchwarz.h
Go to the documentation of this file.
1 /*
2 //@HEADER
3 // ***********************************************************************
4 //
5 // Ifpack: Object-Oriented Algebraic Preconditioner Package
6 // Copyright (2002) Sandia Corporation
7 //
8 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9 // license for use of this work by or on behalf of the U.S. Government.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39 //
40 // ***********************************************************************
41 //@HEADER
42 */
43 
44 #ifndef IFPACK_ADDITIVESCHWARZ_H
45 #define IFPACK_ADDITIVESCHWARZ_H
46 
47 #include "Ifpack_ConfigDefs.h"
48 #include "Ifpack_Preconditioner.h"
49 #include "Ifpack_ConfigDefs.h"
50 #include "Ifpack_Preconditioner.h"
51 #include "Ifpack_Reordering.h"
52 #include "Ifpack_RCMReordering.h"
53 #include "Ifpack_METISReordering.h"
54 #include "Ifpack_LocalFilter.h"
55 #include "Ifpack_NodeFilter.h"
56 #include "Ifpack_SingletonFilter.h"
57 #include "Ifpack_ReorderFilter.h"
58 #include "Ifpack_Utils.h"
60 #include "Epetra_CombineMode.h"
61 #include "Epetra_MultiVector.h"
62 #include "Epetra_Map.h"
63 #include "Epetra_Comm.h"
64 #include "Epetra_Time.h"
65 #include "Epetra_LinearProblem.h"
66 #include "Epetra_RowMatrix.h"
67 #include "Epetra_CrsMatrix.h"
68 #include "Teuchos_ParameterList.hpp"
69 #include "Teuchos_RefCountPtr.hpp"
70 
71 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
72 #include "Ifpack_SubdomainFilter.h"
73 #include "EpetraExt_Reindex_CrsMatrix.h"
74 #include "EpetraExt_Reindex_MultiVector.h"
75 #endif
76 #ifdef IFPACK_NODE_AWARE_CODE
77 #include "EpetraExt_OperatorOut.h"
78 #include "EpetraExt_RowMatrixOut.h"
79 #include "EpetraExt_BlockMapOut.h"
80 #endif
81 
82 #ifdef HAVE_IFPACK_AMESOS
83  #include "Ifpack_AMDReordering.h"
84 #endif
85 
86 
88 
141 template<typename T>
143 
144 public:
145 
147 
160  int OverlapLevel_in = 0);
161 
165 
167 
169 
178  virtual int SetUseTranspose(bool UseTranspose_in);
180 
182 
184 
194  virtual int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
195 
197 
208  virtual int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
209 
211  virtual double NormInf() const;
213 
215 
217  virtual const char * Label() const;
218 
220  virtual bool UseTranspose() const;
221 
223  virtual bool HasNormInf() const;
224 
226  virtual const Epetra_Comm & Comm() const;
227 
229  virtual const Epetra_Map & OperatorDomainMap() const;
230 
232  virtual const Epetra_Map & OperatorRangeMap() const;
234 
236  virtual bool IsInitialized() const
237  {
238  return(IsInitialized_);
239  }
240 
242  virtual bool IsComputed() const
243  {
244  return(IsComputed_);
245  }
246 
248 
268 
269  // @}
270 
271  // @{ Query methods
272 
274  virtual int Initialize();
275 
277  virtual int Compute();
278 
280  virtual double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
281  const int MaxIters = 1550,
282  const double Tol = 1e-9,
283  Epetra_RowMatrix* Matrix_in = 0);
284 
286  virtual double Condest() const
287  {
288  return(Condest_);
289  }
290 
292  virtual const Epetra_RowMatrix& Matrix() const
293  {
294  return(*Matrix_);
295  }
296 
298  virtual bool IsOverlapping() const
299  {
300  return(IsOverlapping_);
301  }
302 
304  virtual std::ostream& Print(std::ostream&) const;
305 
306  virtual const T* Inverse() const
307  {
308  return(&*Inverse_);
309  }
310 
312  virtual int NumInitialize() const
313  {
314  return(NumInitialize_);
315  }
316 
318  virtual int NumCompute() const
319  {
320  return(NumCompute_);
321  }
322 
324  virtual int NumApplyInverse() const
325  {
326  return(NumApplyInverse_);
327  }
328 
330  virtual double InitializeTime() const
331  {
332  return(InitializeTime_);
333  }
334 
336  virtual double ComputeTime() const
337  {
338  return(ComputeTime_);
339  }
340 
342  virtual double ApplyInverseTime() const
343  {
344  return(ApplyInverseTime_);
345  }
346 
348  virtual double InitializeFlops() const
349  {
350  return(InitializeFlops_);
351  }
352 
353  virtual double ComputeFlops() const
354  {
355  return(ComputeFlops_);
356  }
357 
358  virtual double ApplyInverseFlops() const
359  {
360  return(ApplyInverseFlops_);
361  }
362 
364  virtual int OverlapLevel() const
365  {
366  return(OverlapLevel_);
367  }
368 
370  virtual const Teuchos::ParameterList& List() const
371  {
372  return(List_);
373  }
374 
375 protected:
376 
377  // @}
378 
379  // @{ Internal merhods.
380 
383  { }
384 
386  int Setup();
387 
388  // @}
389 
390  // @{ Internal data.
391 
393  Teuchos::RefCountPtr<const Epetra_RowMatrix> Matrix_;
395  Teuchos::RefCountPtr<Ifpack_OverlappingRowMatrix> OverlappingMatrix_;
397 /*
398  //TODO if we choose to expose the node aware code, i.e., no ifdefs,
399  //TODO then we should switch to this definition.
400  Teuchos::RefCountPtr<Epetra_RowMatrix> LocalizedMatrix_;
401 */
402 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
406  Teuchos::RCP<Epetra_CrsMatrix> SubdomainCrsMatrix_;
408  Teuchos::RCP<Epetra_CrsMatrix> ReindexedCrsMatrix_;
409 
410  // The following data members are needed when doing ApplyInverse
411  // with an Ifpack_SubdomainFilter local matrix
412  Teuchos::RCP<Epetra_Map> tempMap_;
413  Teuchos::RCP<Epetra_Map> tempDomainMap_;
414  Teuchos::RCP<Epetra_Map> tempRangeMap_;
415  Teuchos::RCP<EpetraExt::CrsMatrix_Reindex> SubdomainMatrixReindexer_;
416  Teuchos::RCP<EpetraExt::MultiVector_Reindex> DomainVectorReindexer_;
417  Teuchos::RCP<EpetraExt::MultiVector_Reindex> RangeVectorReindexer_;
418  mutable Teuchos::RCP<Epetra_MultiVector> tempX_;
419  mutable Teuchos::RCP<Epetra_MultiVector> tempY_;
420 #else
421 # ifdef IFPACK_NODE_AWARE_CODE
422  Teuchos::RefCountPtr<Ifpack_NodeFilter> LocalizedMatrix_;
423 # else
424  Teuchos::RefCountPtr<Ifpack_LocalFilter> LocalizedMatrix_;
425 # endif
426 #endif
427 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
428  // Some data members used for determining the subdomain id (color)
430  int MpiRank_;
432  int NumMpiProcs_;
434  int NumMpiProcsPerSubdomain_;
436  int NumSubdomains_;
438  int SubdomainId_;
439 #endif
440  std::string Label_;
457  double Condest_;
463  std::string ReorderingType_;
465  Teuchos::RefCountPtr<Ifpack_Reordering> Reordering_;
467  Teuchos::RefCountPtr<Ifpack_ReorderFilter> ReorderedLocalizedMatrix_;
471  Teuchos::RefCountPtr<Ifpack_SingletonFilter> SingletonFilter_;
477  mutable int NumApplyInverse_;
481  double ComputeTime_;
483  mutable double ApplyInverseTime_;
489  mutable double ApplyInverseFlops_;
491  Teuchos::RefCountPtr<Epetra_Time> Time_;
493  Teuchos::RefCountPtr<T> Inverse_;
495 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
496  mutable Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingX;
497  mutable Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingY;
498 #endif
499 # ifdef IFPACK_NODE_AWARE_CODE
500  mutable Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingX;
501  mutable Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingY;
502 #endif
503 
504 }; // class Ifpack_AdditiveSchwarz<T>
505 
506 //==============================================================================
507 template<typename T>
510  int OverlapLevel_in) :
511 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
512  MpiRank_(0),
513  NumMpiProcs_(1),
514  NumMpiProcsPerSubdomain_(1),
515  NumSubdomains_(1),
516  SubdomainId_(0),
517 #endif
518  IsInitialized_(false),
519  IsComputed_(false),
520  UseTranspose_(false),
521  IsOverlapping_(false),
522  OverlapLevel_(OverlapLevel_in),
523  CombineMode_(Zero),
524  Condest_(-1.0),
525  ComputeCondest_(true),
526  UseReordering_(false),
527  ReorderingType_("none"),
528  FilterSingletons_(false),
529  NumInitialize_(0),
530  NumCompute_(0),
531  NumApplyInverse_(0),
532  InitializeTime_(0.0),
533  ComputeTime_(0.0),
534  ApplyInverseTime_(0.0),
535  InitializeFlops_(0.0),
536  ComputeFlops_(0.0),
537  ApplyInverseFlops_(0.0)
538 {
539  // Construct a reference-counted pointer with the input matrix, don't manage the memory.
540  Matrix_ = Teuchos::rcp( Matrix_in, false );
541 
542  if (Matrix_->Comm().NumProc() == 1)
543  OverlapLevel_ = 0;
544 
545  if ((OverlapLevel_ != 0) && (Matrix_->Comm().NumProc() > 1))
546  IsOverlapping_ = true;
547  // Sets parameters to default values
548  Teuchos::ParameterList List_in;
549  SetParameters(List_in);
550 }
551 
552 #ifdef IFPACK_NODE_AWARE_CODE
553 extern int ML_NODE_ID;
554 #endif
555 
556 //==============================================================================
557 template<typename T>
559 {
560  using std::cerr;
561  using std::endl;
562 
563  Epetra_RowMatrix* MatrixPtr;
564 
565 #ifndef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
566 # ifdef IFPACK_NODE_AWARE_CODE
567 /*
568  sleep(3);
569  if (Comm().MyPID() == 0) cout << "Printing out ovArowmap" << endl;
570  Comm().Barrier();
571 
572  EpetraExt::BlockMapToMatrixMarketFile("ovArowmap",OverlappingMatrix_->RowMatrixRowMap());
573  if (Comm().MyPID() == 0) cout << "Printing out ovAcolmap" << endl;
574  Comm().Barrier();
575  EpetraExt::BlockMapToMatrixMarketFile("ovAcolmap",OverlappingMatrix_->RowMatrixColMap());
576  Comm().Barrier();
577 */
578 /*
579  EpetraExt::RowMatrixToMatlabFile("ovA",*OverlappingMatrix_);
580  fprintf(stderr,"p %d n %d matrix file done\n",Comm().MyPID(),ML_NODE_ID);
581  Comm().Barrier();
582 */
583  int nodeID;
584  try{ nodeID = List_.get("ML node id",0);}
585  catch(...){fprintf(stderr,"%s","Ifpack_AdditiveSchwarz<T>::Setup(): no parameter \"ML node id\"\n\n");
586  std::cout << List_ << std::endl;}
587 # endif
588 #endif
589 
590  try{
591  if (OverlappingMatrix_ != Teuchos::null)
592  {
593 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
594  if (NumMpiProcsPerSubdomain_ > 1) {
595  LocalizedMatrix_ = Teuchos::rcp(new Ifpack_SubdomainFilter(OverlappingMatrix_, SubdomainId_));
596  } else {
597  LocalizedMatrix_ = Teuchos::rcp(new Ifpack_LocalFilter(OverlappingMatrix_));
598  }
599 #else
600 # ifdef IFPACK_NODE_AWARE_CODE
601  Ifpack_NodeFilter *tt = new Ifpack_NodeFilter(OverlappingMatrix_,nodeID); //FIXME
602  LocalizedMatrix_ = Teuchos::rcp(tt);
603  //LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(OverlappingMatrix_) );
604 # else
605  LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(OverlappingMatrix_) );
606 # endif
607 #endif
608  }
609  else
610  {
611 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
612 
613  if (NumMpiProcsPerSubdomain_ > 1) {
614  LocalizedMatrix_ = Teuchos::rcp(new Ifpack_SubdomainFilter(Matrix_, SubdomainId_));
615  } else {
616  LocalizedMatrix_ = Teuchos::rcp(new Ifpack_LocalFilter(Matrix_));
617  }
618 #else
619 # ifdef IFPACK_NODE_AWARE_CODE
620  Ifpack_NodeFilter *tt = new Ifpack_NodeFilter(Matrix_,nodeID); //FIXME
621  LocalizedMatrix_ = Teuchos::rcp(tt);
622  //LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(Matrix_) );
623 # else
624  LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(Matrix_) );
625 # endif
626 #endif
627  }
628  }
629  catch(...) {
630  fprintf(stderr,"%s","AdditiveSchwarz Setup: problem creating local filter matrix.\n");
631  }
632 
633  if (LocalizedMatrix_ == Teuchos::null)
634  IFPACK_CHK_ERR(-5);
635 
636  // users may want to skip singleton check
637  if (FilterSingletons_) {
638  SingletonFilter_ = Teuchos::rcp( new Ifpack_SingletonFilter(LocalizedMatrix_) );
639  MatrixPtr = &*SingletonFilter_;
640  }
641  else
642  MatrixPtr = &*LocalizedMatrix_;
643 
644  if (UseReordering_) {
645 
646  // create reordering and compute it
647  if (ReorderingType_ == "rcm")
648  Reordering_ = Teuchos::rcp( new Ifpack_RCMReordering() );
649  else if (ReorderingType_ == "metis")
650  Reordering_ = Teuchos::rcp( new Ifpack_METISReordering() );
651 #ifdef HAVE_IFPACK_AMESOS
652  else if (ReorderingType_ == "amd" )
653  Reordering_ = Teuchos::rcp( new Ifpack_AMDReordering() );
654 #endif
655  else {
656  cerr << "reordering type not correct (" << ReorderingType_ << ")" << endl;
657  exit(EXIT_FAILURE);
658  }
659  if (Reordering_ == Teuchos::null) IFPACK_CHK_ERR(-5);
660 
661  IFPACK_CHK_ERR(Reordering_->SetParameters(List_));
662  IFPACK_CHK_ERR(Reordering_->Compute(*MatrixPtr));
663 
664  // now create reordered localized matrix
665  ReorderedLocalizedMatrix_ =
666  Teuchos::rcp( new Ifpack_ReorderFilter(Teuchos::rcp( MatrixPtr, false ), Reordering_) );
667 
668  if (ReorderedLocalizedMatrix_ == Teuchos::null) IFPACK_CHK_ERR(-5);
669 
670  MatrixPtr = &*ReorderedLocalizedMatrix_;
671  }
672 
673 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
674  // The subdomain matrix needs to be reindexed by Amesos so we need to make a CrsMatrix
675  // and then reindex it with EpetraExt.
676  // The reindexing is done here because this feature is only implemented in Amesos_Klu,
677  // and it is desired to have reindexing with other direct solvers in the Amesos package
678 
679  SubdomainCrsMatrix_.reset(new Epetra_CrsMatrix(Copy, MatrixPtr->RowMatrixRowMap(), -1));
680  Teuchos::RCP<Epetra_Import> tempImporter = Teuchos::rcp(new Epetra_Import(SubdomainCrsMatrix_->Map(), MatrixPtr->Map()));
681  SubdomainCrsMatrix_->Import(*MatrixPtr, *tempImporter, Insert);
682  SubdomainCrsMatrix_->FillComplete();
683 
684  if (NumMpiProcsPerSubdomain_ > 1) {
685  tempMap_.reset(new Epetra_Map(SubdomainCrsMatrix_->RowMap().NumGlobalElements(),
686  SubdomainCrsMatrix_->RowMap().NumMyElements(),
687  0, SubdomainCrsMatrix_->Comm()));
688  tempRangeMap_.reset(new Epetra_Map(SubdomainCrsMatrix_->OperatorRangeMap().NumGlobalElements(),
689  SubdomainCrsMatrix_->OperatorRangeMap().NumMyElements(),
690  0, SubdomainCrsMatrix_->Comm()));
691  tempDomainMap_.reset(new Epetra_Map(SubdomainCrsMatrix_->OperatorDomainMap().NumGlobalElements(),
692  SubdomainCrsMatrix_->OperatorDomainMap().NumMyElements(),
693  0, SubdomainCrsMatrix_->Comm()));
694 
695  SubdomainMatrixReindexer_.reset(new EpetraExt::CrsMatrix_Reindex(*tempMap_));
696  DomainVectorReindexer_.reset(new EpetraExt::MultiVector_Reindex(*tempDomainMap_));
697  RangeVectorReindexer_.reset(new EpetraExt::MultiVector_Reindex(*tempRangeMap_));
698 
699  ReindexedCrsMatrix_.reset(&((*SubdomainMatrixReindexer_)(*SubdomainCrsMatrix_)), false);
700 
701  MatrixPtr = &*ReindexedCrsMatrix_;
702 
703  Inverse_ = Teuchos::rcp( new T(&*ReindexedCrsMatrix_) );
704  } else {
705  Inverse_ = Teuchos::rcp( new T(&*SubdomainCrsMatrix_) );
706  }
707 #else
708  Inverse_ = Teuchos::rcp( new T(MatrixPtr) );
709 #endif
710 
711  if (Inverse_ == Teuchos::null)
712  IFPACK_CHK_ERR(-5);
713 
714  return(0);
715 }
716 
717 //==============================================================================
718 template<typename T>
720 {
721 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
722  MpiRank_ = Matrix_->Comm().MyPID();
723  NumMpiProcs_ = Matrix_->Comm().NumProc();
724  NumMpiProcsPerSubdomain_ = List_in.get("subdomain: number-of-processors", 1);
725  NumSubdomains_ = NumMpiProcs_ / NumMpiProcsPerSubdomain_;
726  SubdomainId_ = MpiRank_ / NumMpiProcsPerSubdomain_;
727 
728  if (NumSubdomains_ == 1) {
729  IsOverlapping_ = false;
730  }
731 
732 #endif
733  // compute the condition number each time Compute() is invoked.
734  ComputeCondest_ = List_in.get("schwarz: compute condest", ComputeCondest_);
735  // combine mode
736  if( Teuchos::ParameterEntry *combineModeEntry = List_in.getEntryPtr("schwarz: combine mode") )
737  {
738  if( typeid(std::string) == combineModeEntry->getAny().type() )
739  {
740  std::string mode = List_in.get("schwarz: combine mode", "Add");
741  if (mode == "Add")
742  CombineMode_ = Add;
743  else if (mode == "Zero")
744  CombineMode_ = Zero;
745  else if (mode == "Insert")
746  CombineMode_ = Insert;
747  else if (mode == "InsertAdd")
748  CombineMode_ = InsertAdd;
749  else if (mode == "Average")
750  CombineMode_ = Average;
751  else if (mode == "AbsMax")
752  CombineMode_ = AbsMax;
753  else
754  {
756  true,std::logic_error
757  ,"Error, The (Epetra) combine mode of \""<<mode<<"\" is not valid! Only the values"
758  " \"Add\", \"Zero\", \"Insert\", \"InsertAdd\", \"Average\", and \"AbsMax\" are accepted!"
759  );
760  }
761  }
762  else if ( typeid(Epetra_CombineMode) == combineModeEntry->getAny().type() )
763  {
764  CombineMode_ = Teuchos::any_cast<Epetra_CombineMode>(combineModeEntry->getAny());
765  }
766  else
767  {
768  // Throw exception with good error message!
769  Teuchos::getParameter<std::string>(List_in,"schwarz: combine mode");
770  }
771  }
772  else
773  {
774  // Make the default be a std::string to be consistent with the valid parameters!
775  List_in.get("schwarz: combine mode","Zero");
776  }
777  // type of reordering
778  ReorderingType_ = List_in.get("schwarz: reordering type", ReorderingType_);
779  if (ReorderingType_ == "none")
780  UseReordering_ = false;
781  else
782  UseReordering_ = true;
783  // if true, filter singletons. NOTE: the filtered matrix can still have
784  // singletons! A simple example: upper triangular matrix, if I remove
785  // the lower node, I still get a matrix with a singleton! However, filter
786  // singletons should help for PDE problems with Dirichlet BCs.
787  FilterSingletons_ = List_in.get("schwarz: filter singletons", FilterSingletons_);
788 
789  // This copy may be needed by Amesos or other preconditioners.
790  List_ = List_in;
791 
792  return(0);
793 }
794 
795 //==============================================================================
796 template<typename T>
798 {
799  IsInitialized_ = false;
800  IsComputed_ = false; // values required
801  Condest_ = -1.0; // zero-out condest
802 
803  if (Time_ == Teuchos::null)
804  Time_ = Teuchos::rcp( new Epetra_Time(Comm()) );
805 
806  Time_->ResetStartTime();
807 
808  // compute the overlapping matrix if necessary
809  if (IsOverlapping_) {
810 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
811  if (NumMpiProcsPerSubdomain_ > 1) {
812  OverlappingMatrix_ = Teuchos::rcp( new Ifpack_OverlappingRowMatrix(Matrix_, OverlapLevel_, SubdomainId_) );
813  } else {
814  OverlappingMatrix_ = Teuchos::rcp( new Ifpack_OverlappingRowMatrix(Matrix_, OverlapLevel_));
815  }
816 #else
817 # ifdef IFPACK_NODE_AWARE_CODE
818  int myNodeID;
819  try{ myNodeID = List_.get("ML node id",-1);}
820  catch(...){fprintf(stderr,"pid %d: no such entry (returned %d)\n",Comm().MyPID(),myNodeID);}
821 /*
822  cout << "pid " << Comm().MyPID()
823  << ": calling Ifpack_OverlappingRowMatrix with myNodeID = "
824  << myNodeID << ", OverlapLevel_ = " << OverlapLevel_ << endl;
825 */
826  OverlappingMatrix_ = Teuchos::rcp( new Ifpack_OverlappingRowMatrix(Matrix_, OverlapLevel_, myNodeID) );
827 # else
828  OverlappingMatrix_ = Teuchos::rcp( new Ifpack_OverlappingRowMatrix(Matrix_, OverlapLevel_) );
829 # endif
830 #endif
831 
832  if (OverlappingMatrix_ == Teuchos::null) {
833  IFPACK_CHK_ERR(-5);
834  }
835  }
836 
837 # ifdef IFPACK_NODE_AWARE_CODE
838 /*
839  sleep(1);
840  Comm().Barrier();
841 */
842 # endif
843 
844  IFPACK_CHK_ERR(Setup());
845 
846 # ifdef IFPACK_NODE_AWARE_CODE
847 /*
848  sleep(1);
849  Comm().Barrier();
850 */
851 #endif
852 
853  if (Inverse_ == Teuchos::null)
854  IFPACK_CHK_ERR(-5);
855 
856  if (LocalizedMatrix_ == Teuchos::null)
857  IFPACK_CHK_ERR(-5);
858 
859  IFPACK_CHK_ERR(Inverse_->SetUseTranspose(UseTranspose()));
860  IFPACK_CHK_ERR(Inverse_->SetParameters(List_));
861  IFPACK_CHK_ERR(Inverse_->Initialize());
862 
863  // Label is for Aztec-like solvers
864  Label_ = "Ifpack_AdditiveSchwarz, ";
865  if (UseTranspose())
866  Label_ += ", transp";
867  Label_ += ", ov = " + Ifpack_toString(OverlapLevel_)
868  + ", local solver = \n\t\t***** `" + std::string(Inverse_->Label()) + "'";
869 
870  IsInitialized_ = true;
871  ++NumInitialize_;
872  InitializeTime_ += Time_->ElapsedTime();
873 
874 #ifdef IFPACK_FLOPCOUNTERS
875  // count flops by summing up all the InitializeFlops() in each
876  // Inverse. Each Inverse() can only give its flops -- it acts on one
877  // process only
878  double partial = Inverse_->InitializeFlops();
879  double total;
880  Comm().SumAll(&partial, &total, 1);
881  InitializeFlops_ += total;
882 #endif
883 
884  return(0);
885 }
886 
887 //==============================================================================
888 template<typename T>
890 {
891  if (IsInitialized() == false)
892  IFPACK_CHK_ERR(Initialize());
893 
894  Time_->ResetStartTime();
895  IsComputed_ = false;
896  Condest_ = -1.0;
897 
898  IFPACK_CHK_ERR(Inverse_->Compute());
899 
900  IsComputed_ = true; // need this here for Condest(Ifpack_Cheap)
901  ++NumCompute_;
902  ComputeTime_ += Time_->ElapsedTime();
903 
904 #ifdef IFPACK_FLOPCOUNTERS
905  // sum up flops
906  double partial = Inverse_->ComputeFlops();
907  double total;
908  Comm().SumAll(&partial, &total, 1);
909  ComputeFlops_ += total;
910 #endif
911 
912  // reset the Label
913  std::string R = "";
914  if (UseReordering_)
915  R = ReorderingType_ + " reord, ";
916 
917  if (ComputeCondest_)
918  Condest(Ifpack_Cheap);
919 
920  // add Condest() to label
921  Label_ = "Ifpack_AdditiveSchwarz, ov = " + Ifpack_toString(OverlapLevel_)
922  + ", local solver = \n\t\t***** `" + std::string(Inverse_->Label()) + "'"
923  + "\n\t\t***** " + R + "Condition number estimate = "
924  + Ifpack_toString(Condest_);
925 
926  return(0);
927 }
928 
929 //==============================================================================
930 template<typename T>
932 {
933  // store the flag -- it will be set in Initialize() if Inverse_ does not
934  // exist.
935  UseTranspose_ = UseTranspose_in;
936 
937  // If Inverse_ exists, pass it right now.
938  if (Inverse_!=Teuchos::null)
939  IFPACK_CHK_ERR(Inverse_->SetUseTranspose(UseTranspose_in));
940  return(0);
941 }
942 
943 //==============================================================================
944 template<typename T>
947 {
948  IFPACK_CHK_ERR(Matrix_->Apply(X,Y));
949  return(0);
950 }
951 
952 //==============================================================================
953 template<typename T>
955 {
956  return(-1.0);
957 }
958 
959 //==============================================================================
960 template<typename T>
962 {
963  return(Label_.c_str());
964 }
965 
966 //==============================================================================
967 template<typename T>
969 {
970  return(UseTranspose_);
971 }
972 
973 //==============================================================================
974 template<typename T>
976 {
977  return(false);
978 }
979 
980 //==============================================================================
981 template<typename T>
983 {
984  return(Matrix_->Comm());
985 }
986 
987 //==============================================================================
988 template<typename T>
990 {
991  return(Matrix_->OperatorDomainMap());
992 }
993 
994 //==============================================================================
995 template<typename T>
997 {
998  return(Matrix_->OperatorRangeMap());
999 }
1000 
1001 //==============================================================================
1002 template<typename T>
1005 {
1006  // compute the preconditioner is not done by the user
1007  if (!IsComputed())
1008  IFPACK_CHK_ERR(-3);
1009 
1010  int NumVectors = X.NumVectors();
1011 
1012  if (NumVectors != Y.NumVectors())
1013  IFPACK_CHK_ERR(-2); // wrong input
1014 
1015  Time_->ResetStartTime();
1016 
1017  Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingX;
1018  Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingY;
1019  Teuchos::RefCountPtr<Epetra_MultiVector> Xtmp;
1020 
1021  // for flop count, see bottom of this function
1022 #ifdef IFPACK_FLOPCOUNTERS
1023  double pre_partial = Inverse_->ApplyInverseFlops();
1024  double pre_total;
1025  Comm().SumAll(&pre_partial, &pre_total, 1);
1026 #endif
1027 
1028  // process overlap, may need to create vectors and import data
1029  if (IsOverlapping()) {
1030 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
1031  if (OverlappingX == Teuchos::null) {
1032  OverlappingX = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(), X.NumVectors()) );
1033  if (OverlappingX == Teuchos::null) IFPACK_CHK_ERR(-5);
1034  } else assert(OverlappingX->NumVectors() == X.NumVectors());
1035  if (OverlappingY == Teuchos::null) {
1036  OverlappingY = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(), Y.NumVectors()) );
1037  if (OverlappingY == Teuchos::null) IFPACK_CHK_ERR(-5);
1038  } else assert(OverlappingY->NumVectors() == Y.NumVectors());
1039 #else
1040 # ifdef IFPACK_NODE_AWARE_CODE
1041  if (OverlappingX == Teuchos::null) {
1042  OverlappingX = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
1043  X.NumVectors()) );
1044  if (OverlappingX == Teuchos::null) IFPACK_CHK_ERR(-5);
1045  } else assert(OverlappingX->NumVectors() == X.NumVectors());
1046  if (OverlappingY == Teuchos::null) {
1047  OverlappingY = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
1048  Y.NumVectors()) );
1049  if (OverlappingY == Teuchos::null) IFPACK_CHK_ERR(-5);
1050  } else assert(OverlappingY->NumVectors() == Y.NumVectors());
1051 #else
1052  OverlappingX = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
1053  X.NumVectors()) );
1054  OverlappingY = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
1055  Y.NumVectors()) );
1056  if (OverlappingY == Teuchos::null) IFPACK_CHK_ERR(-5);
1057 # endif
1058 #endif
1059  OverlappingY->PutScalar(0.0);
1060  OverlappingX->PutScalar(0.0);
1061  IFPACK_CHK_ERR(OverlappingMatrix_->ImportMultiVector(X,*OverlappingX,Insert));
1062  // FIXME: this will not work with overlapping and non-zero starting
1063  // solutions. The same for other cases below.
1064  // IFPACK_CHK_ERR(OverlappingMatrix_->ImportMultiVector(Y,*OverlappingY,Insert));
1065  }
1066  else {
1067  Xtmp = Teuchos::rcp( new Epetra_MultiVector(X) );
1068  OverlappingX = Xtmp;
1069  OverlappingY = Teuchos::rcp( &Y, false );
1070  }
1071 
1072  if (FilterSingletons_) {
1073  // process singleton filter
1074  Epetra_MultiVector ReducedX(SingletonFilter_->Map(),NumVectors);
1075  Epetra_MultiVector ReducedY(SingletonFilter_->Map(),NumVectors);
1076  IFPACK_CHK_ERR(SingletonFilter_->SolveSingletons(*OverlappingX,*OverlappingY));
1077  IFPACK_CHK_ERR(SingletonFilter_->CreateReducedRHS(*OverlappingY,*OverlappingX,ReducedX));
1078 
1079  // process reordering
1080  if (!UseReordering_) {
1081  IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReducedX,ReducedY));
1082  }
1083  else {
1084  Epetra_MultiVector ReorderedX(ReducedX);
1085  Epetra_MultiVector ReorderedY(ReducedY);
1086  IFPACK_CHK_ERR(Reordering_->P(ReducedX,ReorderedX));
1087  IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReorderedX,ReorderedY));
1088  IFPACK_CHK_ERR(Reordering_->Pinv(ReorderedY,ReducedY));
1089  }
1090 
1091  // finish up with singletons
1092  IFPACK_CHK_ERR(SingletonFilter_->UpdateLHS(ReducedY,*OverlappingY));
1093  }
1094  else {
1095  // process reordering
1096  if (!UseReordering_) {
1097 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
1098  if (NumMpiProcsPerSubdomain_ > 1) {
1099  tempX_.reset(&((*RangeVectorReindexer_)(*OverlappingX)), false);
1100  tempY_.reset(&((*DomainVectorReindexer_)(*OverlappingY)), false);
1101  IFPACK_CHK_ERR(Inverse_->ApplyInverse(*tempX_,*tempY_));
1102  } else {
1103  IFPACK_CHK_ERR(Inverse_->ApplyInverse(*OverlappingX, *OverlappingY));
1104  }
1105 #else
1106  IFPACK_CHK_ERR(Inverse_->ApplyInverse(*OverlappingX,*OverlappingY));
1107 #endif
1108  }
1109  else {
1110  Epetra_MultiVector ReorderedX(*OverlappingX);
1111  Epetra_MultiVector ReorderedY(*OverlappingY);
1112  IFPACK_CHK_ERR(Reordering_->P(*OverlappingX,ReorderedX));
1113  IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReorderedX,ReorderedY));
1114  IFPACK_CHK_ERR(Reordering_->Pinv(ReorderedY,*OverlappingY));
1115  }
1116  }
1117 
1118  if (IsOverlapping()) {
1119  IFPACK_CHK_ERR(OverlappingMatrix_->ExportMultiVector(*OverlappingY,Y,
1120  CombineMode_));
1121  }
1122 
1123 #ifdef IFPACK_FLOPCOUNTERS
1124  // add flops. Note the we only have to add the newly counted
1125  // flops -- and each Inverse returns the cumulative sum
1126  double partial = Inverse_->ApplyInverseFlops();
1127  double total;
1128  Comm().SumAll(&partial, &total, 1);
1129  ApplyInverseFlops_ += total - pre_total;
1130 #endif
1131 
1132  // FIXME: right now I am skipping the overlap and singletons
1133  ++NumApplyInverse_;
1134  ApplyInverseTime_ += Time_->ElapsedTime();
1135 
1136  return(0);
1137 
1138 }
1139 
1140 //==============================================================================
1141 template<typename T>
1142 std::ostream& Ifpack_AdditiveSchwarz<T>::
1143 Print(std::ostream& os) const
1144 {
1145  using std::endl;
1146 
1147 #ifdef IFPACK_FLOPCOUNTERS
1148  double IF = InitializeFlops();
1149  double CF = ComputeFlops();
1150  double AF = ApplyInverseFlops();
1151 
1152  double IFT = 0.0, CFT = 0.0, AFT = 0.0;
1153  if (InitializeTime() != 0.0)
1154  IFT = IF / InitializeTime();
1155  if (ComputeTime() != 0.0)
1156  CFT = CF / ComputeTime();
1157  if (ApplyInverseTime() != 0.0)
1158  AFT = AF / ApplyInverseTime();
1159 #endif
1160 
1161  if (Matrix().Comm().MyPID())
1162  return(os);
1163 
1164  os << endl;
1165  os << "================================================================================" << endl;
1166  os << "Ifpack_AdditiveSchwarz, overlap level = " << OverlapLevel_ << endl;
1167  if (CombineMode_ == Insert)
1168  os << "Combine mode = Insert" << endl;
1169  else if (CombineMode_ == Add)
1170  os << "Combine mode = Add" << endl;
1171  else if (CombineMode_ == Zero)
1172  os << "Combine mode = Zero" << endl;
1173  else if (CombineMode_ == Average)
1174  os << "Combine mode = Average" << endl;
1175  else if (CombineMode_ == AbsMax)
1176  os << "Combine mode = AbsMax" << endl;
1177 
1178  os << "Condition number estimate = " << Condest_ << endl;
1179  os << "Global number of rows = " << Matrix_->NumGlobalRows64() << endl;
1180 
1181 #ifdef HAVE_IFPACK_PARALLEL_SUBDOMAIN_SOLVERS
1182  os << endl;
1183  os << "================================================================================" << endl;
1184  os << "Subcommunicator stats" << endl;
1185  os << "Number of MPI processes in simulation: " << NumMpiProcs_ << endl;
1186  os << "Number of subdomains: " << NumSubdomains_ << endl;
1187  os << "Number of MPI processes per subdomain: " << NumMpiProcsPerSubdomain_ << endl;
1188 #endif
1189 
1190  os << endl;
1191  os << "Phase # calls Total Time (s) Total MFlops MFlops/s" << endl;
1192  os << "----- ------- -------------- ------------ --------" << endl;
1193  os << "Initialize() " << std::setw(5) << NumInitialize()
1194  << " " << std::setw(15) << InitializeTime()
1195 #ifdef IFPACK_FLOPCOUNTERS
1196  << " " << std::setw(15) << 1.0e-6 * IF
1197  << " " << std::setw(15) << 1.0e-6 * IFT
1198 #endif
1199  << endl;
1200  os << "Compute() " << std::setw(5) << NumCompute()
1201  << " " << std::setw(15) << ComputeTime()
1202 #ifdef IFPACK_FLOPCOUNTERS
1203  << " " << std::setw(15) << 1.0e-6 * CF
1204  << " " << std::setw(15) << 1.0e-6 * CFT
1205 #endif
1206  << endl;
1207  os << "ApplyInverse() " << std::setw(5) << NumApplyInverse()
1208  << " " << std::setw(15) << ApplyInverseTime()
1209 #ifdef IFPACK_FLOPCOUNTERS
1210  << " " << std::setw(15) << 1.0e-6 * AF
1211  << " " << std::setw(15) << 1.0e-6 * AFT
1212 #endif
1213  << endl;
1214  os << "================================================================================" << endl;
1215  os << endl;
1216 
1217  return(os);
1218 }
1219 
1220 #include "Ifpack_Condest.h"
1221 //==============================================================================
1222 template<typename T>
1224 Condest(const Ifpack_CondestType CT, const int MaxIters,
1225  const double Tol, Epetra_RowMatrix* Matrix_in)
1226 {
1227  if (!IsComputed()) // cannot compute right now
1228  return(-1.0);
1229 
1230  Condest_ = Ifpack_Condest(*this, CT, MaxIters, Tol, Matrix_in);
1231 
1232  return(Condest_);
1233 }
1234 
1235 #endif // IFPACK_ADDITIVESCHWARZ_H
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the preconditioner to X, returns the result in Y.
int Setup()
Sets up the localized matrix and the singleton filter.
int OverlapLevel_
Level of overlap among the processors.
std::string ReorderingType_
Type of reordering of the local matrix.
Teuchos::RefCountPtr< const Epetra_RowMatrix > Matrix_
Pointers to the matrix to be preconditioned.
virtual const Epetra_Map & RowMatrixRowMap() const =0
double ComputeFlops_
Contains the number of flops for Compute().
virtual int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied (not implemented).
virtual const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
Ifpack_AMDReordering: approximate minimum degree reordering.
T & get(ParameterList &l, const std::string &name)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
bool UseTranspose_
If true, solve with the transpose (not supported by all solvers).
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
bool ComputeCondest_
If true, compute the condition number estimate each time Compute() is called.
const int NumVectors
Definition: performance.cpp:71
double Condest_
Contains the estimated condition number.
virtual bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
Ifpack_METISReordering: A class to reorder a graph using METIS.
double ApplyInverseFlops_
Contain sthe number of flops for ApplyInverse().
virtual bool IsComputed() const
Returns true if the preconditioner has been successfully computed.
Teuchos::RefCountPtr< Ifpack_OverlappingRowMatrix > OverlappingMatrix_
Pointers to the overlapping matrix.
virtual int OverlapLevel() const
Returns the level of overlap.
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
virtual bool IsOverlapping() const
Returns true is an overlapping matrix is present.
virtual double ComputeTime() const
Returns the time spent in Compute().
Teuchos::RefCountPtr< Ifpack_Reordering > Reordering_
Pointer to a reordering object.
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
Teuchos::RefCountPtr< Ifpack_ReorderFilter > ReorderedLocalizedMatrix_
Pointer to the reorderd matrix.
ParameterEntry * getEntryPtr(const std::string &name)
Ifpack_CondestType
Ifpack_CondestType: enum to define the type of condition number estimate.
std::string Label_
Contains the label of this object.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
bool IsOverlapping_
If true, overlapping is used.
virtual std::ostream & Print(std::ostream &) const
Prints major information about this preconditioner.
virtual int Initialize()
Initialized the preconditioner.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Ifpack_AdditiveSchwarz(const Ifpack_AdditiveSchwarz &RHS)
Copy constructor (should never be used)
double ComputeTime_
Contains the time for all successful calls to Compute().
Ifpack_RCMReordering: reverse Cuthill-McKee reordering.
Ifpack_AdditiveSchwarz: a class to define Additive Schwarz preconditioners of Epetra_RowMatrix&#39;s.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
virtual const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
virtual const Epetra_BlockMap & Map() const =0
virtual bool UseTranspose() const
Returns the current UseTranspose setting.
Ifpack_AdditiveSchwarz(Epetra_RowMatrix *Matrix_in, int OverlapLevel_in=0)
Ifpack_AdditiveSchwarz constructor with given Epetra_RowMatrix.
double InitializeFlops_
Contains the number of flops for Initialize().
virtual bool HasNormInf() const
Returns true if the this object can provide an approximate Inf-norm, false otherwise.
Ifpack_SingletonFilter: Filter based on matrix entries.
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the matrix to X, returns the result in Y.
virtual ~Ifpack_AdditiveSchwarz()
Destructor.
Ifpack_OverlappingRowMatrix: matrix with ghost rows, based on Epetra_RowMatrix.
Ifpack_LocalFilter a class for light-weight extraction of the submatrix corresponding to local rows a...
bool IsInitialized_
If true, the preconditioner has been successfully initialized.
std::string Ifpack_toString(const int &x)
Converts an integer to std::string.
#define false
int NumInitialize_
Contains the number of successful calls to Initialize().
virtual int SetParameters(Teuchos::ParameterList &List)
Sets the parameters.
bool UseReordering_
If true, reorder the local matrix.
double Ifpack_Condest(const Ifpack_Preconditioner &IFP, const Ifpack_CondestType CT, const int MaxIters, const double Tol, Epetra_RowMatrix *Matrix)
virtual int Compute()
Computes the preconditioner.
int NumApplyInverse_
Contains the number of successful call to ApplyInverse().
Teuchos::RefCountPtr< Ifpack_SingletonFilter > SingletonFilter_
filtering object.
double InitializeTime_
Contains the time for all successful calls to Initialize().
cheap estimate
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
bool IsComputed_
If true, the preconditioner has been successfully computed.
virtual const Teuchos::ParameterList & List() const
Returns a reference to the internally stored list.
Epetra_CombineMode
virtual const char * Label() const
Returns a character string describing the operator.
Teuchos::RefCountPtr< T > Inverse_
Pointer to the local solver.
int NumCompute_
Contains the number of successful call to Compute().
virtual const T * Inverse() const
virtual double InitializeTime() const
Returns the time spent in Initialize().
virtual const Epetra_Comm & Comm() const
Returns a pointer to the Epetra_Comm communicator associated with this operator.
double ApplyInverseTime_
Contains the time for all successful calls to ApplyInverse().
#define IFPACK_CHK_ERR(ifpack_err)
Ifpack_ReorderFilter: a class for light-weight reorder of local rows and columns of an Epetra_RowMatr...
virtual double NormInf() const
Returns the infinity norm of the global matrix (not implemented)
Teuchos::ParameterList List_
Stores a copy of the list given in SetParameters()
#define true
Epetra_CombineMode CombineMode_
Combine mode for off-process elements (only if overlap is used)
virtual const Epetra_RowMatrix & Matrix() const
Returns a refernence to the internally stored matrix.
bool FilterSingletons_
Filter for singletons.
virtual double Condest() const
Returns the estimated condition number, or -1.0 if not computed.
virtual int NumCompute() const
Returns the number of calls to Compute().
Teuchos::RefCountPtr< Epetra_Time > Time_
Object used for timing purposes.
#define RHS(a)
Definition: MatGenFD.c:60
Teuchos::RefCountPtr< Ifpack_LocalFilter > LocalizedMatrix_
Localized version of Matrix_ or OverlappingMatrix_.