1 / 25

ML: Multilevel Preconditioning Package

ML: Multilevel Preconditioning Package. Trilinos User’s Group Meeting Wednesday, October 15, 2003 Jonathan Hu. 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. Outline.

norman
Download Presentation

ML: Multilevel Preconditioning Package

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. ML: Multilevel Preconditioning Package Trilinos User’s Group Meeting Wednesday, October 15, 2003 Jonathan Hu 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 • What is ML? • Basic multigrid concepts • Configuring and building ML • Interoperability with other packages • User options available within ML • Example: Using ML as a preconditioner • What to do if something goes wrong • Sandia apps that use ML • Future Plans • Documentation, mailing lists, and getting help

  3. ML Package • C package that provides multigrid preconditioning for linear solver methods • Developers: Ray Tuminaro, Jonathan Hu, Charles Tong (LLNL) • Main methods • Geometric • Grid refinement hierarchy • 2-level FE basis function domain decomposition • AMG (algebraic) • Smoothed aggregation • Edge-element AMG for Maxwell’s equations • Classical AMG

  4. Solve A4u4 = f4 (S4) I3 R3 Approximate PDE on (user supplied) grid hierarchy A3u3 = f3 (S3) R2 I2 A2u2 = f2 (S2) R1 I1 A1u1 = f1 (S1) Develop smoothers Sk (approximate solve on a level) Jacobi, Gauss-Seidel, CG, etc. Develop grid transfer operators: restriction Rk and interpolation Ik. Geometric Multigrid Use coarse solves (k<4) to accelerate convergence for A4

  5. Algebraic Multigrid (AMG) Same Goal: Solve ALuL=fL. • Build MG operators (Ak, Ik, Rk, Sk’s) automatically to define a hierarchy: Akuk=fk, k=1,…,L • Main difference: In AMG, interpolation operators Ik’s built automatically • also main difficulty • Once Ik’s are defined, the rest follows “easily”: • Rk = IkT (usually) • Ak-1 = Rk-1 Ak Ik-1 (triple matrix product) • Smoother (iterative method) Sk • Gauss-Seidel, polynomial, conjugate gradient, etc.

  6. W Cycle FMV Cycle Recursive MG Algorithm (V Cycle) V Cycle to solve A4u4=f4 MG(f, u, k) { if (k == 1) u1 = (A1)-1f1 else { Sk(Ak, fk,uk) //pre-smooth rk = fk – Ak uk fk-1 = Rk-1 rk; uk-1 = 0 //restrict MG(fk-1, uk-1, k-1) uk = uk + Ik-1 uk-1 //interpolate & correct Sk(Ak, fk,uk) //post-smooth } } Ak, Rk, Ik, Sk required on all levels Sk: smoothers Rk: restriction operators Ik: interpolation operators k=4 k=1

  7. Algebraic construction of Ik : Aggregation • Greedy algorithm • parallel can be complicated

  8. 5 5 2 2 7 2 7 5 7 2 { 1 if ith point within jth aggregate 0 otherwise Algebraic Construction of Ik: Coefficients • Build tentativetIkto interpolate null space where tIk(i,j) = Finding Ik • Smoothed aggregation • Improves tIk with Jacobi’s method: Ik = (I -  diag(Ak) -1 Ak) tIk • Ik emphasizes what is not smoothed by Jacobi

  9. ML Capabilities • MG cycling: V, W, full V, full W • Grid Transfers • Several automatic coarse grid generators • Several automatic grid transfer operators • Smoothers • Jacobi, Gauss-Seidel, Hiptmair, LU, Aztec methods, sparse approximate inverses, polynomial • Kernels: matrix/matrix multiply, etc.

  10. Configuring and Building ML • ML builds by default when Trilinos is configured/built • Configure help • See ml/README file, or ml/configure --help • Enabling direct solver support (default is off) configure --with-ml_superlu \ --with-incdirs=“-I/usr/local/superlu/include” \ --with-ldflags=“/usr/local/superlu/libsuperlu.a” • Performance monitoring --enable-ml_timing --enable_ml_flops • Example suite builds by default Your-build-location/packages/ml/examples

  11. Via Epetra & TSF Accepts user data as Epetra objects Aztecoo TSF Can be wrapped as Epetra_Operator ML Meros TSF interface exists Epetra Amesos Amesos interface coming very soon Other solvers Accepts other solvers and MatVecs Other matvecs ML Interoperability with Trilinos Packages

  12. Common Decisions for ML Users • What smoother to use • # of smoothing steps • Coarsening strategy • Cycling strategy

  13. ML Smoother Choices • Jacobi • Simplest, cheapest, usually least effective. • Damping parameter () needed • Point Gauss-Seidel • Equation satisfied one unknown at a time • Better than Jacobi • May need damping • Can be problematic in parallel • processor-based (stale off-proc values)

  14. ML Smoother Choices • Block Gauss-Seidel • Satisfy several equations simultaneously by modifying several DOFs (inverting subblock associated with DOFs). • Blocks can correspond to aggregates • Aztec smoothers • Any Aztec preconditioner can be used as a smoother. • Probably most interesting is ILU & ILUT methods: may be more robust than Gauss-Seidel. • Sparse LU • Hiptmair • Specialized 2-stage smoother for Maxwell’s Eqns. • MLS • Approximation to inverse based on Chebyshev polynomials of smoothed operator. • Competitive with true Gauss-Seidel in serial. • Doesn’t degrade with # processors (unlike processor-based GS)

  15. ML Cycling Choices • V is default (usually works best). • W more expensive, may be more robust. • Full MG (V cycle) more expensive • Less conventional within preconditioners These choices decide how frequently coarse grids are visited compared to fine grid.

  16. ML Aggregation Choices • MIS • Expensive, usually best with many processors • Uncoupled (UC) • Cheap, usually works well. • Hybrid = UC + MIS • Uncoupled on fine grids, MIS on coarser grids

  17. A Simple Example: Trilinos/Packages/aztecoo/example/MLAztecOO Solve Ax=b: A: linear elasticity problem Linear solver: Conjugate Gradient Precond.: AMG smoothed aggregation Solution component Example methods Packages used Linear Solver CG AztecOO (Epetra) ML: multi- grid pre- cond. AMG ML

  18. Simple Example #include “ml_include.h” #include “ml_epetra_operator.h” #include “ml_epetra_utils.h” Epetra_CrsMatrix A; Epetra_Vector x,bb; ... //matrix, vectors loaded here // Construct Epetra Linear Problem Epetra_LinearProblem problem(&A, &x, &bb); // Construct a solver object AztecOO solver(problem); solver.SetAztecOption(AZ_solver, AZ_cg); // Create and set an ML multilevel preconditioner int N_levels = 10; // max # of multigrid levels possible ML_Set_PrintLevel(3); // how much ML info is output to screen ML *ml_handle; ML_Create(&ml_handle,N_levels); // Make linear operator A accessible to ML EpetraMatrix2MLMatrix(ml_handle, 0, &A);

  19. Simple Example (contd.) // Create multigrid hierarchy ML_Aggregate *agg_object; ML_Aggregate_Create(&agg_object); ML_Aggregate_Set_CoarsenScheme_Uncoupled(agg_object); N_levels = ML_Gen_MGHierarchy_UsingAggregation(ml_handle, 0, ML_INCREASING, agg_object); // Set symmetric Gauss-Seidel smoother for MG method int nits = 1; double dampingfactor = ML_DDEFAULT; ML_Gen_Smoother_SymGaussSeidel(ml_handle, ML_ALL_LEVELS, ML_BOTH, nits, dampingfactor); ML_Gen_Solver(ml_handle, ML_MGV, 0, N_levels-1); // Set preconditioner within Epetra Epetra_ML_Operator MLop(ml_handle,comm,map,map); solver.SetPrecOperator(&MLop); // Set some Aztec solver options and iterate solver.SetAztecParam(AZ_rthresh, 1.4); solver.SetAztecParam(AZ_athresh, 10.0); int Niters = 500; solver.Iterate(Niters, 1e-12);

  20. Multigrid: Issues to be aware of • Severely stretched grids or anisotropies • Loss of diagonal dominance • Atypical stencils • Jumps in material properties • Non-symmetric matrices • Boundary conditions • Systems of PDEs • Non-trivial null space (Maxwell’s equations, elasticity)

  21. Things to Try if Multigrid Isn’t Working • Smoothers • Vary number of smoothing steps. • More `robust’ smoothers: block Gauss-Seidel, ILU, MLS polynomial (especially if degradation in parallel). • Vary damping parameters (smaller is more conservative). • Try different aggregation schemes • Try fewer levels • Try drop tolerances, if … • high complexity (printed out by ML). • Severely stretched grids • Anisotropic problems • Variable regions • Reduce prolongator damping parameter • Concerned about operators properties (e.g. highly non-symmetric).

  22. ALEGRA/NEVADA Zpinch Results CG preconditioned with V(1,1) AMG cycle ||r||2/||r0||2 < 10-8 2-stage Hiptmair smoother (4th order polys in each stage) Edge element interpolation

  23. Flow and Transport Solution Algorithms for Chemical Attack in an Airport Terminal: 3D prototype Algorithm; Steady State: Newton – Krylov (GMRES)

  24. Future Improvements • Amesos interface • More direct solvers • Self-correcting capabilities • Better analysis tools • Improved coarse grids • “On-the-fly” adjustments • “Re-partitioning”

  25. Getting Help • See Trilinos/packages/ml/examples • See guide: • ML User’s Guide, ver. 2.0 (in ml/doc) • ML User’s Guide, ver. 3.0 (under construction) • Mailing lists • ml-users@software.sandia.gov • ml-announce@software.sandia.gov • Bug reporting, enhancement requests via bugzilla: • http://software.sandia.gov/bugzilla • Email us directly • jhu@sandia.gov • tuminaro@ca.sandia.gov

More Related