An Overview of Meros

1. meros An Overview of Meros Trilinos User’s Group Wednesday, November 2, 2005 Victoria Howle Computational Sciences and Mathematics Research Department (8962) 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 Meros? • Motivation & background • Incompressible Navier–Stokes • Block preconditioners • Some preconditioners being developed in Meros • A few results from these methods • Code example: user level • Code example: inside Meros • Release plans, etc. • References

3. What is Meros? • Segregated preconditioner package in Trilinos • Scalable block preconditioning for problems that couple simultaneous solution variables • Initial focus is on (incompressible) Navier-Stokes • Release version in progress • Updating (from old TSF) to Thyra interface • Plan to release next Fall ‘06 • Team • Ray Tuminaro 1414, Computational Mathematics & Algorithms • Robert Shuttleworth Univ. of Maryland, Summer Student Intern 2003, 2004, 2005 • Other collaborators • Howard Elman, University of Maryland • Jacob Schroder, University of Illinois, Summer Intern 2005 • John Shadid, Sandia, NM • David Silvester, Manchester Univerity

4. Analyst or Designer Optimization Code (e.g., APPSPACK, MOOCHO, Opt++, Split) Simulation Code (e.g., MPSalsa, Sierra, Sundance) Nonlinear Solver (e.g., NOX) Linear Solver (e.g., AZTECOO, Belos) Preconditioner (e.g., IFPACK, Meros, ML) Linear Algebra Kernel (e.g., Epetra) Where is Meros in The Big Picture • The speed, scalability, and robustness of an application can be heavily dependent on the speed, scalability, and robustness of the linear solvers • Linear algebra often accounts for >80% of the computational time in many applications • Iterative linear solvers are essential in ASC-scale problems • Preconditioning is the key to iterative solver performance

5. Incompressible Navier–Stokes • Examples of incompressible flow problems • Airflow in an airport; e.g., transport of an airborne toxin • Chemical Vapor Deposition • Goal: efficient and robust solution of steady and transient chemically reacting flow applications • Current testbed application: MPSalsa • Early user: Sundance • Related Sandia applications: • Charon • ARIA • Fuego • … Airport source detection problem CVD Reactor

6. Incompressible Navier–Stokes •  = 0 ) steady state,  = 1 ) transient • (2,2)-block = 0 (unstabilized) or = C (stabilized) • Incompressibility constraint ) difficult for linear solvers • Chemically reactive flow  multiphysics; even harder • Indefinite, strongly coupled, nonlinear, nonsymmetric systems

7. Block preconditioners • Want the scalability of multigrid (mesh-independence) • Difficult to apply multigrid to the whole system • Solution: • Segregate blocks and apply multigrid separately to subproblems • Consider the following class of preconditioners: • is an optimal (right) preconditioner when is the Schur complement, (Assuming C = 0)

8. Choosing (Kay & Loghin, Fp) • Key is choosing a good Schur complement approximation to • Motivation: move F-1 so that it does not appear between GT and G • Suppose we have an Fp such that Then And Giving(Kay, Loghin, & Wathen; Silvester, Elman, Kay & Wathen) (Ap is pressure Poisson)

9. Other Choices for • Kay & Loghin Fp method works well, but… • Fp is not a standard operator for apps(pressure convection–diffusion) • Can be difficult for many applications to provide • Even if they can provide it, they don’t really want to • Other options for : Algebraic pressure convection–diffusion methods: • Sparse Approximate Commutator (SPAC) • Least Squares Commutator (LSC) • Algebraically determine an operator Fp such that

10. Algebraic Commutators • Build Fp column by column via ideas similar to sparse approximate inverses (e.g., Grote & Huckle) ) Sparse Approximate Commutators (SPAC) • Fp is no longer a pressure convection-diffusion operator • Minimize via normal equations ) Least Squares Commutators (LSC) • The tilde’s are hiding an issue of algebraic vs. differential commuting

11. Stabilized LSC (C ≠ 0) • Certain discretizations require stabilization • Stabilization term C • For certain discretizations (GTG) is unstable • Blows up on high frequencies • C built to stabilize • Preconditioner also needs stabilization • In 3 places • Use C for preconditioner stabilization, too

12. Fp vs. DD results:Flow over a diamond in MPSalsa • Linear solve timings • Steady state (harder than transient for linear algebra) • Parallel (on Sandia’s ICC) • Re = 25 • Using development version of Meros hooked into MPSalsa through NOX

13. (Matlab) Results: Fp, LSC, and SPAC • Linear iterations for backward facing step problem on underlying 64x192 grid, Q2-Q1 (stable) discretization. • Linear iterations for backward facing step problem on underlying 128x384 grid, Q2-Q1 (stable) discretization. • Results from ifiss • Academic software package that incorporates our new methods and a few other methods (Elman, Silvester, Ramage)

14. (Matlab) Results: Fp, LSC, Stabilized LSC • Linear iterations for lid driven cavity problem on 32x32 grid, Q1-Q1 (needs stabilization) discretization. • Results from ifiss

15. Transition promising academic methods into methods for ASC applications • Promising methods have been developed • We have extended these methods mathematically to suit more realistic needs • Removing need for nonstandard operators • Stabilization • Currently, software for these methods is mostly in academic (Matlab) codes • Now need to develop software to make them available to more real-world apps through Trilinos

16. Meros • Initial focus is on preconditioners for Navier-Stokes • A number of solvers are being incorporated: • Pressure convection-diffusion preconditioners (today’s focus) • Fp (Kay & Loghin) • LSC (and stabilized LSC) • (SPAC?) • Pressure-projection methodsE.g., SIMPLE (SIMPLEC, SIMPLER, etc.)

17. Trilinos packages in an MPSalsa example Time Loop Component Methods Package Nonlinear Loop MPSalsa Finite Element Epetra Linear Solver Nonlinear Solver Newton-Krylov Methods NOX Block Precond Linear Solver GMRESR Aztec00 (Epetra, Thyra) block preconditioner Meros (Thyra) End NonLin Loop End Time Loop F-1 : GMRES/AMG Ŝ-1 : CG/AMG Aztec00, ML Epetra

18. Meros & Trilinos • Meros is a package within Trilinos • Meros is also a user of many other Trilinos packages • Depends on: • Thyra • Teuchos • (Epetra) • Currently uses: • AztecOO • IFPACK • ML • Could use: • Belos • Amesos • …

19. Example preconditioner:First set up abstract solvers for inner solves // WARNING: Assuming TSF-style handles and assuming I have typedeffed // to hide the Templating // WARNING: Examples include functionality that is not yet available in Thyra // Meros builds a PreconditionerFactory so we can pass it to an abstract linear solver // E.g., K & L preconditioner needs the saddlepoint matrix A, plus Fp and Ap, // and choices of solvers for F and Ap // Inner F solver options: Teuchos::ParameterList FParams; FParams.set(“Solver”, “GMRES”); FParams.set(“Preconditioner”, “ML”); FParams.set(“Max Iters”, 200); FParams.set(“Tolerance”, 1.0e-8); // etc LinearSolver FSolver = new AztecSolver(FParams); // Inner Ap solver options: ApParams.set(“Solver”, “PCG”); ApParams.set(“Preconditioner”, “ML”); // etc LinearSolver ApSolver = new AztecSolver(ApParams);

20. Next set up Schur complement approx. and build the preconditioner // Set up Schur complement approx factory (with solvers if necessary) SchurFactory sfac = new KayLoghinSchurFactory(ApSolver); // Build preconditioner factory with these choices PreconditionerFactory pfac = new KayLoghinFactory(outerMaxIters, outerTol, FSolver, sfac, …) // Group operators that are needed by preconditioner OperatorSource opSrc = new KayLoghinOperatorSource(saddleA, Fp, Ap); // Use preconditioner factory directly in an abstract solver outerParams.set(“Solver”,”GMRESR”); // etc LinearSolver solver = new AztecSolver(outerParams); SolverState solverstate = solver.solve(pfac, opSrc, rhs, soln);

21. Example (cont.) // Get Thyra Preconditioner from factory for a particular set of ops Preconditioner Pinv = pfac.createPreconditioner(opSrc); // Get Thyra LinearOpWithSolve to use precond op more directly LinearOperator Minv = Pinv.right(); outerParams.set(“Solver”,”GMRESR”); LinearSolver solver = new AztecSolver(outerParams); SolverState solverstate = solver.solve(A*Minv, rhs, intermediateSoln); soln = Minv * intermediateSoln; // Simple constructors will make intelligent choices of defaults: PreconditionerFactory pfac = new KayLoghinFactory(maxIters, Tol); // Still need the appropriate operators for the chosen method // (some can be built algebraically by default if not given, e.g., SPAC) OperatorSource opSrc = new OperatorSource(A, Fp, Ap);

22. Inside createPreconditioner() // Build the preconditioner given 2x2 block matrix (etc.) Preconditioner KayLoghinFactory::createPreconditioner(…) { // Get F, G, GT blocks from the block operators LinearOperator F = A.getBlock(1,1); LinearOperator G = A.getBlock(1,2); LinearOperator Gt = A.getBlock(2,1); // LinearOperators Ap and Fp built here or gotten from OpSrc // Set up F solve (given solver and parameters or build with defaults) LinearOpWithSolve Finv = F.inverse(FSolver);

23. createPreconditioner() (cont.) // Setup Schur complement approximation and solver // (given by user or build using defaults) LinearOpWithSolve Apinv = Ap.inverse(ApSolver); LinearOperator Sinv = -Fp * Apinv; // Or if we were building an LSC preconditioner LinearOperator GtG = Gt * G; LinearOpWithSolve GtGinv = Ap.inverse(ApSolver); LinearOperator Sinv = -GtGinv * Gt * F * G * GtGinv;

24. createPreconditioner() (cont.) LinearOperator Iv = IdentityOperator(F.domain()); // velocity space LinearOperator Ip = IdentityOperator(G.domain()); // pressure space // Domain and range of A are Thyra product spaces, velocity x pressure LinearOperator P1 = new BlockLinearOp(A.domain(),A.range()); LinearOperator P2 = new BlockLinearOp(A.domain(),A.range()); LinearOperator P3 = new BlockLinearOp(A.domain(),A.range()); P1.setBlock(1,1,Finv); P1.setBlock(2,2,Ip); P2.setBlock(1,1,Iv); P2.setBlock(2,2,Ip); P2.setBlock(1,2,G); P3.setBlock(1,1,Iv); P3.setBlock(2,2,Sinv); return new GenericRightPreconditioner(P1*P2*P3); }

25. Plans & Info • Planning to release Meros 1.0 in Fall ’06 (with closest major Trilinos release) • Initial block preconditioner selection should include: • Pressure Convection-Diffusion • Kay & Loghin (Fp) • Least Squares Commutator (LSC) • SPAC? • Pressure Projection • SIMPLE • SIMPLEC, SIMPLER? • Web page: software.sandia.gov/Trilinos/packages/meros/index.html • Mailing lists: Meros-Announce, Meros-Users, etc. • vehowle@sandia.gov

26. References • Elman, Silvester, and Wathen, Performance and analysis of saddle point preconditioners for the discrete steady-state Navier-Stokes equations, Numer. Math., 90 (2002), pp. 665-688. • Kay, Loghin, and Wathen, A preconditioner for the steady-state Navier-Stokes equations, SIAM J. Sci. Comput., 2002. • Elman, H., Shadid, and Tuminaro, A Parallel Block Multi-level Preconditioner for the 3D Incompressible Navier-Stokes Equations, J. Comput. Phys, Vol. 187, pp. 504-523, May 2003. • Elman, H., Shadid, Shuttleworth, and Tuminaro, Block Preconditioners Based on Approximate Commutators, to appear in SIAM J. Sci. Comput., Copper Mountain Special Issue, 2005. • Elman, H., Shadid, and Tuminaro, Least Squares Preconditioners for Stabilized Discretizations of the Navier-Stokes Equations, in progress.