1 / 28

Outline

Belos: A Framework for Next-generation Iterative Linear Solvers April 7 th , 2008 Mike Heroux Rich Lehoucq Mike Parks Heidi Thornquist (Lead).

lalasa
Download Presentation

Outline

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Belos: A Framework for Next-generation Iterative Linear SolversApril 7th, 2008Mike HerouxRich LehoucqMike ParksHeidi Thornquist (Lead) Sandia is a multiprogram laboratory operated by Sandia Corporation, a Lockheed Martin Company,for the United States Department of Energy under contract DE-AC04-94AL85000.

  2. Outline • Belos Introduction • Motivation • Belos design • What’s in Trilinos 8.0 • Framework implementations • Current linear solver iterations • Parameter list driven solver managers • Concluding Remarks • What’s next for Belos … • High-level Linear Solver Strategies • Stratimikos

  3. Introducing Belos • Next-generation linear solver library, written in templated C++. • Provide a generic framework for developing algorithms for solving large-scale linear problems. • Algorithm implementation is accomplished through the use of traits classes and abstract base classes: • ex.: MultiVecTraits, OperatorTraits • ex.: Iteration, SolverManager, StatusTest • Includes block linear solvers: • Higher operator performance.

  4. Why Belos? • AztecOO provides solvers for Ax=b • What about solvers for: • Simultaneously solved systems w/ multiple-RHS: AX = B • Sequentially solved systems w/ multiple-RHS: AXi = Bi , i=1,…,t • Sequences of multiple-RHS systems: AiXi = Bi , i=1,…,t • Many advanced methods for these types of linear systems • Block methods: block GMRES [Vital], block CG/BICG [O’Leary] • “Seed” solvers: hybrid GMRES [Nachtigal, et al.] • Recycling solvers: recycled Krylov methods [Parks, et al.] • Restarting techniques, orthogonalization techiques, … • Efficient R&D of advanced linear solution methods requires: • Interoperability • Extensibility • Reusability

  5. Belos Design Design Goal: Create a linear solver library that: • Provides an abstract interface to an operator-vector products, scaling, and preconditioning. • Allows the user to enlist any linear algebra package for the elementary vector space operations essential to the algorithm. (Epetra, PETSc, Thyra, etc.) • Allows the user to define convergence of any algorithm (a.k.a. status testing). • Allows the user to determine the verbosity level, formatting, and processor for the output. • Allows for easier R&D of new solvers through “managers” using “iterations” as the basic kernels.

  6. Status Test Linear Problem Generic Operator Interface Orthogonalization Manager Generic MultiVector Interface Output Manager Generic Iteration Kernel void iterationA.iterate() { … while ( statusTest.checkStatus(this)!= Passed ) { … % Compute operator-vector product linProb->applyOp( P, AP ); % Compute  = PTAP MVT::MvTransMv( 1.0, P, AP,  ); … % Orthonormalize new direction vectors int rank = orthoMgr->normalize( *P ); … outputMgr->print(Belos::Debug, “something”); } }

  7. Solution Strategy Generic Iteration Kernel Generic Solver Manager Belos::ReturnType solverMgrA.solve() { … <initialize iteration> … void iterationA.iterate() { … while ( statusTest.checkStatus(this)!= Passed ) { … % Compute operator-vector product linProb->applyOp( P, AP ); % Compute  = PTAP MVT::MvTransMv( 1.0, P, AP,  ); … % Orthonormalize new direction vectors int rank = orthoMgr->normalize( *P ); … outputMgr->print(Belos::Debug, “something”); } } … <restart, deflate converged solutions> … }

  8. Belos Classes • MultiVecTraits and OperatorTraits • Traits classes for interfacing linear algebra • LinearProblem Class • Describes the problem and stores the answer • OrthoManagerClass • Provide basic interface for orthogonalization • StatusTest Class • Control testing of convergence, etc. • OutputManager Class • Control verbosity and printing in a MP scenario • Iteration Class • Provide basic linear solver iteration interface. • SolverManagerClass • Parameter list driven strategy object describing behavior of solver

  9. What’s in Trilinos 8.0

  10. Linear Algebra Interface • Belos::MultiVecTraits<ST,MV> • Interface to define the linear algebra required by most iterative linear solvers: • creational methods • dot products, norms • update methods • initialize / randomize • Implementations: • MultiVecTraits<double,Epetra_MultiVector> • MultiVecTraits<ST,Thyra::MultiVectorBase<ST> > • Belos::OperatorTraits<ST,MV,OP> • Interface to enable the application of an operator to a multivector. • Implementations: • OperatorTraits<double,Epetra_MultiVector,Epetra_Operator> • OperatorTraits<ST,Thyra::MultiVectorBase<ST>,Thyra::LinearOpBase<ST> >

  11. Linear Problem Interface • Provides an interface between the basic iterations and the linear problem to be solved. • Templated classBelos::LinearProblem<ST,MV,OP> • Allows preconditioning to be removed from the algorithm. • Behavior defined through traits mechanisms. • Methods: • setOperator(…) / getOperator() • setLHS(…) / getLHS() • setRHS(…) / getRHS() • setLeftPrec(…) / getLeftPrec() / isLeftPrec() • setRightPrec(…) / getRightPrec() / isRightPrec() • apply(…) / applyOp(…) / applyLeftPrec(…) / applyRightPrec(…) • setHermitian(…) / isHermitian() • setProblem(…)

  12. Orthogonalization Manager • Abstract interface to orthogonalization / orthonormalization routines for multivectors. • Abstract base class Belos::[Mat]OrthoManager<ST,MV,OP> • void innerProd(…) const; • void norm(…) const; • void project(…) const; • int normalize(…) const; • int projectAndNormalize(…) const; • magnitudeType orthogError(…) const; • magnitudeType orthonormError(…) const; • Implementations: • Belos::DGKSOrthoManager • Belos::ICGSOrthoManager • Belos::IMGSOrthoManager

  13. StatusTest Interface • Informs linear solver iteration of change in state, as defined by user. • Similar to NOX / AztecOO. • Multiple criterion can be logically connected. • Abstract base class Belos::StatusTest<ST,MV,OP> • TestStatus checkStatus( Iteration<…>* iterate ); • TestStatus getStatus() const; • void clearStatus(); • void reset(); • ostream& print( ostream& os, int indent = 0 ) const; • Implementations: • Belos::StatusTestMaxIters • Belos::StatusTestGenResNorm • Belos::StatusTestImpResNorm • Belos::StatusTestOutput • Belos::StatusTestCombo

  14. Output Manager Interface • Templated class that enables control of the linear solver output. • Behavior defined through traits mechanism • Belos::OutputManager<ST> • Get/Set Methods: • void setVerbosity( int vbLevel ); • int getVerbosity( ); • ostream& stream( MsgType type ); • Query Methods: • bool isVerbosity( MsgType type ); • Print Methods: • void print( MsgType type, const string output ); • Message Types: • Belos::Errors, Belos::Warnings, Belos::IterationDetails, Belos::OrthoDetails, Belos::FinalSummary, Belos::TimingDetails, Belos::StatusTestDetails, Belos::Debug • Default is lowest verbosity (Belos::Errors), output on one processor.

  15. Iteration Interface • Provides an abstract interface to Belos basic iteration kernels. • Abstract base class Belos::Iteration<ST,MV,OP> • int getNumIters() const; void resetNumIters(int iter); • Teuchos::RCP<const MV> getNativeResiduals( … ) const; • Teuchos::RCP<const MV> getCurrentUpdate() const; • int getBlockSize() const; void setBlockSize(int blockSize); • const LinearProblem<ST,MV,OP>& getProblem() const; • void iterate(); void initialize(); • Iterations require these classes: • Belos::LinearProblem • Belos::OutputManager • Belos::StatusTest • Belos::OrthoManager • Implementations: • Belos::BlockGmresIter • Belos::BlockFGmresIter • Belos::PseudoBlockGmresIter • Belos::GCRODRIter • Belos::CGIter, • Belos::BlockCGIter

  16. Solver Manager • Provides an abstract interface to Belos solver managers • solver strategies • Abstract base class Belos::SolverManager • void setProblem(…); void setParameters(…); • const Belos::LinearProblem<ST,MV,OP>& getProblem() const; • Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const; • Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters() const; • Belos::ReturnType solve(); • Solvers are parameter list driven, take two arguments: • Belos::LinearProblem • Teuchos::ParameterList [validated] • Implementations: • Belos::BlockGmresSolMgr • Belos::PseudoBlockGmresSolMgr • Belos::GCRODRSolMgr • Belos::BlockCGSolMgr

  17. Examples

  18. #include "BelosConfigDefs.hpp" #include "BelosLinearProblem.hpp" #include "BelosEpetraAdapter.hpp" #include "BelosBlockCGSolMgr.hpp“ … int main(int argc, char *argv[]) { typedef Epetra_MultiVector MV; typedef Epetra_Operator OP; … // Create AX=B Teuchos::RCP<Epetra_Operator> A; Teuchos::RCP<Epetra_MultiVector> X; Teuchos::RCP<Epetra_MultiVector> B; … // Create parameter list Teuchos::ParameterList belosList; belosList.set( "Block Size", blocksize ); belosList.set( "Maximum Iterations", maxiters ); belosList.set( "Convergence Tolerance", tol ); if (verbose) { belosList.set( "Verbosity", Belos::Errors + Belos::Warnings + Belos::TimingDetails + Belos::FinalSummary + Belos::StatusTestDetails ); belosList.set( "Output Frequency", frequency ); } else belosList.set( "Verbosity", Belos::Errors + Belos::Warnings ); // Construct an unpreconditioned linear problem Belos::LinearProblem<double,MV,OP> problem( A, X, B ); problem.setHermitian(); bool set = problem.setProblem(); // VERY IMPORTANT! if (set == false) { return -1; } // Create an iterative solver manager Belos::BlockCGSolMgr<double,MV,OP> newSolver( rcp(&problem,false), rcp(&belosList,false) ); // Perform solve Belos::ReturnType ret = newSolver.solve(); if ( ret!=Belos::Converged ) { // Belos did not converge! return -1; } Block CG Example(belos/example/BlockCG/BlockCGEpetraExFile.cpp)

  19. #include "BelosConfigDefs.hpp" #include "BelosLinearProblem.hpp" #include "BelosEpetraAdapter.hpp" #include "BelosBlockCGSolMgr.hpp“ … int main(int argc, char *argv[]) { typedef Epetra_MultiVector MV; typedef Epetra_Operator OP; … // Create AX=B Teuchos::RCP<OP> A; Teuchos::RCP<MV> X; Teuchos::RCP<MV> B; … // Create parameter list Teuchos::ParameterList belosList; belosList.set( "Block Size", blocksize ); belosList.set( "Maximum Iterations", maxiters ); belosList.set( "Convergence Tolerance", tol ); if (verbose) { belosList.set( "Verbosity", Belos::Errors + Belos::Warnings + Belos::TimingDetails + Belos::FinalSummary + Belos::StatusTestDetails ); belosList.set( "Output Frequency", frequency ); } else belosList.set( "Verbosity", Belos::Errors + Belos::Warnings ); // Construct an SPD preconditioner M // Anything that supports OP interface can go here! Teuchos::RCP<OP> M; … // Construct a preconditioned linear problem Belos::LinearProblem<double,MV,OP> problem( A, X, B ); problem.setHermitian(); if (leftprec) problem.setLeftPrec( M ); else problem.setRightPrec( M ); bool set = problem.setProblem(); // VERY IMPORTANT! if (set == false) { return -1; } // Create an iterative solver manager Belos::BlockCGSolMgr<double,MV,OP> newSolver( rcp(&problem,false), rcp(&belosList,false) ); // Perform solve Belos::ReturnType ret = newSolver.solve(); if ( ret!=Belos::Converged ) { // Belos did not converge! return -1; } Preconditioned Block CG Example(belos/example/BlockCG/BlockPrecCGEpetraExFile.cpp)

  20. #include "BelosConfigDefs.hpp" #include "BelosLinearProblem.hpp" #include "BelosEpetraAdapter.hpp" #include "BelosBlockGmresSolMgr.hpp" … int main(int argc, char *argv[]) { typedef Epetra_MultiVector MV; typedef Epetra_Operator OP; … // Create AX=B Teuchos::RCP<OP> A; Teuchos::RCP<MV> X; Teuchos::RCP<MV> B; … // Create parameter list Teuchos::ParameterList belosList; belosList.set( "Num Blocks", maxsubspace ); belosList.set( "Block Size", blocksize ); belosList.set( "Maximum Iterations", maxiters ); belosList.set( "Maximum Restarts", maxrestarts ); belosList.set( "Convergence Tolerance", tol ); if (verbose) { belosList.set( "Verbosity", Belos::Errors + Belos::Warnings + Belos::TimingDetails + Belos::FinalSummary + Belos::StatusTestDetails ); belosList.set( "Output Frequency", frequency ); } else belosList.set( "Verbosity", Belos::Errors + Belos::Warnings ); // Construct a preconditioner M // Anything that supports OP interface can go here! Teuchos::RCP<OP> M; … // Construct a preconditioned linear problem Belos::LinearProblem<double,MV,OP> problem( A, X, B ); if (leftprec) problem.setLeftPrec( M ); else problem.setRightPrec( M ); bool set = problem.setProblem(); // VERY IMPORTANT! if (set == false) { return -1; } // Create an iterative solver manager Belos::BlockGmresSolMgr<double,MV,OP> newSolver( rcp(&problem,false), rcp(&belosList,false) ); // Perform solve Belos::ReturnType ret = newSolver.solve(); if ( ret!=Belos::Converged ) { // Belos did not converge! if ( newSolver.isLOADetected() ) { // Convergence was impeded by a loss of accuracy } return -1; } Preconditioned Block GMRES Example(belos/example/BlockGmres/BlockPrecGmresEpetraExFile.cpp)

  21. What’s next for Belos? • Belos provides a powerful framework for developing robust linear solvers, but there’s more work to do … • Future Development: • More performance improvements • More advanced solution methods • Check out the Trilinos Tutorial: http://trilinos.sandia.gov/Trilinos8.0Tutorial.pdf • See Belos website for more info: http://trilinos.sandia.gov/packages/belos

  22. High-level Linear SolverStrategies UsingStratimikos

  23. Nonlinear Algorithms and Applications : Everyone for Themselves? Nonlinear ANA Solvers in Trilinos NOX / LOCA … Rythmos MOOCHO Trilinos and non-Trilinos Preconditioner and Linear Solver Capability Sandia Applications … Xyce Charon Tramonto Aria Aleph Key Point • BAD

  24. Introducing Stratimikos • Stratimikos created Greek words "stratigiki“ (strategy) and "grammikos“ (linear) • Based on the foundation of abstract interface layer Thyra • Defines class Stratimikos::DefaultLinearSolverBuilder • Provides common access to: • Linear Solvers: Amesos, AztecOO, Belos, … • Preconditioners: Ifpack, ML, … • Reads in options through a parameter list (read from XML?) • Accepts any linear system objects that provide • Epetra_Operator / Epetra_RowMatrix view of the matrix • SPMD vector views for the RHS and LHS (e.g. Epetra_[Multi]Vector objects) • Provides uniform access to linear solver options that can be leveraged across multiple applications and algorithms • Future: TOPS-2 will add PETSc and other linear solvers and preconditioners! Key Points • Stratimikos is an important building block for creating more sophisticated linear solver capabilities!

  25. Preconditioners and Preconditioner Factories PreconditionerFactoryBase : Creates and initializes PrecondtionerBase objects <<create>> prec PreconditionerFactoryBase createPrec() : PreconditionerBase initializePrec( in fwdOp, inout prec ) PreconditionerBase getLeftPrecOp() : LinearOpBase getRightPrecOp() : LinearOpBase getUnspecifiedPrecOp() : LinearOpBase Create preconditioner prec with preconditioner operators PL and/or PR such that PLA, or APR, or PLAPR is “easier” to solve than unpreconditioned A. • Allows unlimited creation/reuse of preconditioner objects • Supports reuse of factorization structures • Adapters currently available for Ifpack and ML • New Stratimikos package provides a single parameter-driven wrapper for all of these Key Points • You can create your own PreconditionerFactory subclass!

  26. Linear Operator With Solve and Factories LinearOpWithSolveBase : Combines a linear operator and a linear solver LinearOpBase LinearOpWithSolveBase solve( in B, inout X, … ) • Appropriate for both direct and iterative solvers • Supports multiple simultaneous solutions as multi-vectors • Allows targeting of different solution criteria to different RHSs • Supports a “default” solve Key Points • You can create your own subclass! LinearOpWithSolveFactoryBase : Uses LinearOpBase objects to initialize LOWSB objects LinearOpWithSolveFactoryBase createOp() : LinearOpWithSolveBase initializeOp( in fwdOp, inout Op ) initializePreconditionedOp( in fwdOp, in prec, inout Op) <<create>> LinearOpWithSolveBase • Allows unlimited creation/reuse of LinearOpWithSolveBase objects • Supports reuse of factorizations/preconditioners • Supports client-created external preconditioners (which are ignored by direct solvers) • Appropriate for both direct and iterative solvers • Concrete adaptors for Amesos, AztecOO, and Belos are available • New Stratimikos package provides a single parameter-driven wrapper to all of these!

  27. Stratimikos Parameter List and Sublists <ParameterList name=“Stratimikos”> <Parameter name="Linear Solver Type" type="string" value=“AztecOO"/> <Parameter name="Preconditioner Type" type="string" value="Ifpack"/> <ParameterList name="Linear Solver Types"> <ParameterList name="Amesos"> <Parameter name="Solver Type" type="string" value="Klu"/> <ParameterList name="Amesos Settings"> <Parameter name="MatrixProperty" type="string" value="general"/> ... <ParameterList name="Mumps"> ... </ParameterList> <ParameterList name="Superludist"> ... </ParameterList> </ParameterList> </ParameterList> <ParameterList name="AztecOO"> <ParameterList name="Forward Solve"> <Parameter name="Max Iterations" type="int" value="400"/> <Parameter name="Tolerance" type="double" value="1e-06"/> <ParameterList name="AztecOO Settings"> <Parameter name="Aztec Solver" type="string" value="GMRES"/> ... </ParameterList> </ParameterList> ... </ParameterList> <ParameterList name="Belos"> ... </ParameterList> </ParameterList> <ParameterList name="Preconditioner Types"> <ParameterList name="Ifpack"> <Parameter name="Prec Type" type="string" value="ILU"/> <Parameter name="Overlap" type="int" value="0"/> <ParameterList name="Ifpack Settings"> <Parameter name="fact: level-of-fill" type="int" value="0"/> ... </ParameterList> </ParameterList> <ParameterList name="ML"> ... </ParameterList> </ParameterList> </ParameterList> Top level parameters Sublists passed on to package code! Linear Solvers Every parameter and sublist not in red is handled by Thyra code and is fully validated! Preconditioners

  28. Future Plans for Stratimikos • TOPS-2 PETSc/Trilinos Interoperability • Fill PETSc matrices & vectors and use Trilinos preconditioners and linear solvers • Access PETSc preconditioners and solvers through a parameter list in Stratimikos? • Fortran 2003 interface through Stratimikos • Built on Fotran/Epetra wrappers • GUI to manipulate Parameter List XML files • “Smart” linear solver & preconditioner builders?

More Related