Direct and iterative sparse linear solvers applied to groundwater flow simulations

1. Matrix Analysis and Applications October 2007 Direct and iterative sparse linear solvers applied to groundwater flow simulations

2. Jocelyne Erhel INRIA Rennes Jean-Raynald de Dreuzy CNRS, Geosciences Rennes Anthony Beaudoin LMPG, Le Havre Partly funded by Grid’5000 french project

3. Physical context: groundwater flow Tracer evolution during one year (Made, Mississippi) Heterogeneous permeability Dispersion Injection of tracer Flow From Barlebo et al. (2004) Flow governed by the heterogeneous permeability Solute transport by advection and dispersion

4. Numerical modelling strategy Head natural system Characterization of heterogeneity Model validation Simulation of flow and solute transport Numerical Stochastic models Physical model Simulation results

5. Uncertainty Quantification methods Porous geological media Spatial heterogeneity Lack of observations Stochastic models of flow and solute transport -random velocity field -random solute transfer time and dispersivity Horizontal velocity Solute dispersion Permeability field

6. HYDROLAB: parallel software for hydrogeoloy PARALLEL-BASED SCIENTIFIC PLATFORM HYDROLAB Physical models Porous Media Fracture Networks Fractured-Porous Media Numerical methods Utilitaries Input / Output Visualization Results structures Parameters structures Parallel and grid tools Geometry Solvers PDE solvers ODE solvers Linear solvers Particle tracker UQ methods Monte-Carlo Open source libraries Boost, FFTW, CGal, MPI, Hypre, Sundials, OpenGL, Xerces-C Object-oriented and modular with C++ Parallel algorithms with MPI Efficient numerical libraries

7. Physical equations Flow equations: Darcy law and mass conservation Transport equations: advection and dispersion • Saturated medium: one water phase • Constant density: no saltwater • Constant porosity and constant viscosity • Linear equations • Steady-state flow or transient flow • Inert transport: no coupling with chemistry • No coupling between flow and transport • No coupling with heat equations • No coupling with mechanical equations • Classical boundary conditions • Classical initial conditions

8. Flow and transport equations • Flow equations Nul flux and  C/ n=0 • Advection-dispersion equations • Boundary conditions • Initial condition injection Fixed head and  C/ n=0 Fixed head and C=0 Nul flux and  C/ n = 0

9. Monte-Carlo simulations For j=1,…,M generate permeability field Kj using a regular mesh Compute Vj using a finite volume method End For

10. Discrete flow numerical model Linear system Ax=b b: boundary conditions and source term A is a sparse matrix : NZ coefficients Matrix-Vector product : O(NZ) opérations Regular 2D mesh : N=n2 and NZ=5N Regular 3D mesh : N= n3 and NZ=7N Need for parallel sparse linear solvers

11. Accuracy: condition number and variance Estimation with MUMPS solver Cond(A) in O(exp(2)) But in theory, cond(A) in O(Kmax/Kmin) thus in O(exp())

12. Accuracy: condition number and scaling Estimation with Matlab without scaling and with scaling Scaled condition number in O(exp()) As expected

13. Componentwise condition number n=64 Componentwise condition number estimated by || |A-1| |A| |x| + |A-1| b ||1 / ||x||1 Solution error means ||x-xs|| / ||xs||

14. Accuracy: condition number and system size Estimation with MUMPS for =1 Cond(A) in O(N) as expected Condition number not too large for ·3 and for N up to 16 millions

15. Sparse direct linear solver UMFPACK multifrontal solver Robust to variance  but CPU time in O(N1.5) As expected

16. Preconditioned Conjugate Gradient PCG with IC(0) slightly sensitive to variance  But very sensitive to size N Need for a multilevel preconditioner

17. Geometric multigrid HYPRE Solver SMG Linear CPU time in O(N) but sensitivity to variance As expected

18. Algebraic multigrid HYPRE Solver AMG Robust to variance and linear CPU time in O(N) As expected Less efficient than SMG for small variance

19. Algebraic multigrid with 3D domains Robust to variance and CPU time in O(N) Same properties as in 2D

20. Parallel computing facilities Numerical model Clusters at Inria Rennes Grid’5000 project 67.1 millions of unknowns in 3 minutes with 32 processors

21. Parallel performances with 2D domains Parallel CPU time in O(N) SMG more efficient than AMG for small  AMG much more efficient than SMG for large 

22. Macro-dispersion analysis Transversal dispersion Longitudinal dispersion Each curve represents 100 simulations on domains with 67.1 millions of unknowns

23. Conclusion Summary • Efficient and accurate algebraic multigrid solver for groundwater flow in heterogeneous porous media • Good performances with clusters • Macro-dispersion analysis in 2D domains Current and Future work • 3D heterogeneous porous media • Subdomain method with Aitken-Schwarz acceleration • Transient flow in 2D and 3D porous media • Grid computing and parametric simulations • UQ methods