Adaptive Multigrid for Lattice QCD

1. Adaptive Multigrid for Lattice QCD James Brannick October 22, 2007 Penn State brannick@psu.edu

2. Presentation plan • Intro to Lattice Quantum Chromodynamics • Review of Multigrid Basics • (Adaptive) MG for QCD • Numerical Experiments

3. Participants in the MG-QCD project D. Keyes (Columbia University) O. Livne (Univ. of Utah) I. Livshits (Ball State) S. MacLachlan (TU Delft) T. Manteuffel (CU Boulder) S. McCormick (CU Boulder) V. Nistor (PSU) K. Osterlee (TU Delft) J. Osborne (ANL) C. Rebbi (BU) J. Ruge (CU Boulder) P. Varanas (LLNL) A. Bessen (Columbia University) A. Brandt (UCLA, WIS) J. Brannick (PSU) M. Brezina (CU Boulder) R. Brower (BU) M. Clark (BU) R. Falgout (CASC-LLNL) A. Frommer (Wuppertal) K. Kahl (Wuppertal) C. Ketelson (CU Boulder)

4. Forces in Standard Model: SU(N) Nuclei: Weak N=2 (isospin) Atoms: Maxwell N=1(charge) Sub nuclear: Strong N=3 (color) Standard Model: U(1) £ SU(2) £ SU(3)

6. quark QCD path integral Anti-quark Gauge Dirac Operator Generalized Curl (Maxwell)

7. The Dirac PDE (for Quarks) 3x3 color gauge matrices x = (x1 ,x2 ,x3 ,x4 ) (space,time) 4 x 4 sparse spin matrices: 4 non-zero entries 1,-1, i, -i Wilson (1974): discretization on a hypercubic lattice U(x,x+) = exp[i hg A(x)]

8. Discrete Dirac operator on hypercubic lattice Spin projection Operator Color a,b = 1,2,3 Dimension: =1,2,…,d x x+ Spin i,j = 1,2,3,4 x2 axis  x1 axis 

9. Typical lattice size : 323 х 24 • Dimension of M then • 163х 24 х 2 х 3 х 4≈ 106х 106!! • Need M-1(U),Tr[M-1(U)], Det[M(U)] • Together account for dominant cost • more than 80% of the overall flops! • Different types of algorithms used for • different fermionic actions: Krylov • methods typical, as M* = M M(U) • Goal of our NSF PetaApps project is inversion on 2564 lattice - overall • simulation at this resolution requires sustained Petaflop years! Wilson fermion matrix: M /

10. 2-d “toy” problem: Schwinger Model • Space time is 2-d • Gauge links: U(x,x+) = exp[i e A(x)], U(1) theory • Dirac fields have 2 spins (not 4) • Operator is quaternionic (Pauli) matrix involving1, 2

11. Spectrum of Schwinger matrix & “critical slowing down’’ mass gap = .1 mass gap = .01 mass gap = .001

12. Stationary Linear Iterations • Consider solving the linear system using a SLI: • Letek=u - ukbe the error, and note thatrk=Aek. Theerror • propagation operator is then , withrk=f - Auk

13. smoothing Fine Grid Smaller Coarse Grid Multigrid V-cycle Multigrid solvers are optimal (O(N) operations), and hence have good scaling potential prolongation (interpolation) MG uses a sequence of coarse-grid problems to accelerate the solution of the original problem restriction

14. Why does standard MG work? • Low mode is constant -- key to MG success is smooth error, e.g., Laplacian • Constant exactly preserved on coarse level • All near zero modes also preserved!

15. Imaginary Part Real Part Lattice QCD MG • Gauge field and hence low modes not geometrically smooth (locally oscillatory) Lowest mode of M, ¯ = 6 Geometric MG completely fails, preserving low modes for gauge fields requires adaptivity

16. Lattice QCD & (A)MG Lattice QCD and MG have a long and painful history (15+ yrs.) * PTMG (Lauwers et al, 93) * Renormilazation MG (de Forcrand et al, 91) * Projective MG (Brower et al, 91) * Many others .. All failed for non-smooth fields in the m ! mcr limit, failed because did not consider why MG works in first place! Diamond: Jacobi Circle: CG Square: V-cycle Star: W-cycle

17. Adaptive Smoothed Aggregation (SA) Multigrid, Brezina, Falgout, MacLachlan, Manteuffel, McCormick, and Ruge, SIAM Review, 2006. Adaptive Smoothed Aggregation Multigrid for Lattice QCD, Brannick, Brezina, Keyes, Livne, Livshits, MacLachlan, Manteuffel, McCormick, and Ruge, Zikatanov,J. Comp. Physics, 2006. Adaptive Smooth Aggregation AMG Adaptively use (s)low modes to define the evolving V-cycle Initial Algorithm Setup 1. Relax on Ax = 0, with random initial guess (resulting error is repesentative of (s)low modes). 2. Cut vector into pieces over aggregates (blocks). 3. Define the prolongator so that x = Px . 4. Smooth the vector using simple Richardson iteration: P = (I - ¸ A)P, ¸ choosen to minimize condition number of coarse scale operator: A = ((I-¸ A)P)A(I-¸ A)P. c * c General setup adaptive process repeated with current solver to find additional vectors as needed

18. Adaptive Smooth Aggregation MG d s c • Adaptive SA designed for problem without underlying geometry • Uses algebraic defintion of “strength of connection” to define aggregates • QCD defined on regular lattice with unitary connections and so use regular geometric blocking strategy (i.e., 4 x n x n ) • Maintains simple regular geometry on coarse scales, allowing for perfect load balancing and minimal comm. • Original algorithm requires HPD operator and thus we solve M M *

19. * Results for M M • Schwinger model, 2-d with U(1) background on128 x 128 lattice with ¯ = 6,10, Q = 0,4, mass gap = 0.001 - 0.5 • Use 4 x 4 ( x 2 ) blocking and 3 levels with 8 vectors • Under-relaxed MinRes smoother (Bank and Douglas) • Compare MG-PCG with CG

20. Results… ¯ = 6 ¯ = 10 • Critical slowing down eliminated! • No dependence on ¯ • Dramatic improvement over CG Adaptive Multigrid for QCD, Brannick, Brower, Clark, Osborne, and Rebbi, Phys. Rev. Letters, sub. 2007. Adaptive Multigrid for Wilson Fermions, Brannick, Brower, Clark, Osborne, and Rebbi, Proc. of Science: Lattice 07, sub. 2007.

21. Results: MG for M ¯ = 6, N = 128 • Huge reduction of flops • Fewer vectors needed • Develop. of code for full 4-d QCD underway Adaptive Multigrid for the non-hermitian Wilson operator, Brannick, Brower, Clark, Osborne, and Rebbi, in preparation.

22. Possible future collaborations • Current projects (with J. Xu & L. Zikatanov) • Radiation transport • Electromagnetics • Oil resevoir simulations • Fuel cell dynamics • Stochastic matrices • Lattice field theories