Deflated Conjugate Gradient Method for modelling Groundwater Flow in a Layered Grid

1. Deflated Conjugate Gradient Method for modelling Groundwater Flow in a Layered Grid Lennart Ros Deltares & TU Delft Utrecht July 3 2008: 10.30 www.deltares.com Thesis Committee: Prof. Dr. Ir. C. Vuik (TU Delft) Dr M. Genseberger (Deltares) Ir. J. Verkaik (Deltares) Dr H.M. Schuttelaars (TU Delft)

2. Deflated CG method Outline

3. Deflated CG method Outline • Introduction • Deltares • Subsurface, Geohydrology, Clay & Faults • MODFLOW • IBRAHYM & problem • Equation, Discretization & Method • A Simple Example (Matlab) • Deflation Techniques • Results for IBRAHYM • Conclusions & Recommendations • Questions

4. Deflated CG method Introduction

5. Introduction Deltares January1st 2008 Deflated CG method

6. Introduction Subsurface • Subsurface is schematized in layers . • Successive sand and clay (aquifers and aquitards) • Assumption: • Horizontal flow in aquifer • Vertical flow in aquitard Deflated CG method

7. Deflated CG method Introduction Geohydrology • The driving force for groundwater flowis the difference in height and pressure. • To represent this difference we introduce the concept of hydraulicheads, h [L].

8. Deflated CG method Introduction Clay • Very high resistance. • Main property: • Extreme low permeability (vertical) Medium Faults • Vertical barriers inside aquifers. • Main property: • Extreme low permeability (horizontal).

9. Deflated CG method Introduction Clay • Very high resistance. • Main property: • Extreme low permeability (vertical) Consequence: Large contrast in medium parameters Medium Faults • Vertical barriers inside aquifers. • Main property: • Extreme low permeability (horizontal).

10. Deflated CG method Introduction MODFLOW: • MODFLOW is a software package which calculates hydraulic heads. • Developed by theU.S. Geological Survey. • Open-source code: everyone can use and improvethis program • Rectangular grid and uses cell-centered variables. • Quasi-3D model.

11. Deflated CG method Introduction IBRAHYM: • groundwater model developed for several waterboards inLimburg. • uses at most 19 layers to model groundwater flow area. • uses grid cells of 25 times 25 meter to get detailed information. • a lot of clay and faults • these cause model to suffer from bad convergence behavior ofsolver.

12. Deflated CG method Equation, Discretization & Method

13. Deflated CG method Equation, Discretization & Method Governing Equation: Where:

14. Deflated CG method Equation, Discretization & Method Finite Volume Discretization:

15. Deflated CG method Equation, Discretization & Method Discretized Equation Using Finite Volume Method: Where:

16. Deflated CG method Equation, Discretization & Method Solution Method in MODFLOW: stop criteria inner loop: or: maximal number of inner iteration is reached • MODFLOW uses outer and inner iteration loops • We look at inner iteration loop: • solves a linear system of equations: • preconditioner: • iterative method: • MODFLOW uses a: • Modified Incomplete Cholesky Conjugate Gradient Method • with: SOR

17. Deflated CG method A Simple Example

18. Deflated CG method A Simple Example Simple Testcase: • 2 Dimensional Problem • 15 rows, 15 colums • A contrast in the parameter • on 1/3th of the domain

19. Deflated CG method A Simple Example Observations for a simple testcase in Matlab: Assume A has eigenvalues: Preconditioning MODFLOW: Modified Incomplete Cholesky Preconditoning generally works

20. Deflated CG method A Simple Example Observations for a simple testcase in Matlab:

21. Deflated CG method A Simple Example Observations for a simple testcase in Matlab: Smallest eigenvalue: 0.00010283296716 Next eigenvalue: 0.04870854847951

22. Deflated CG method A Simple Example Observations for a simple testcase in Matlab: • Due to the small eigenvalue we have a slow converging model. • Want to get rid of this eigenvalue(s) GOAL: IMPROVE CONVERGENCE BEHAVIOUR OF THE IBRAHYM MODEL IDEA: USE DEFLATION BASED PRECONDITIONER

23. Deflated CG method Deflation Techniques

24. Deflated CG method Deflation Techniques Basic Idea of Deflation: General linear system of equations: Now define: And instead we solve the deflated system:

25. Deflated CG method Deflation Techniques Deflation using Eigenvectors: Assume that A has eigenvalues: and we choose the corresponding eigenvectors such that If we now define Then: PROBLEM: eigenvalues and eigenvectors are generally unknown

26. Deflated CG method Deflation Techniques Alternative Deflation Techniques: • Random Subdomain Deflation • Deflation based on Physics: • Use layers as boundary of domain (1 domain is 1 layer) Original domain Subdomains

27. Deflation Techniques Results for the simple problem: • Deflation using subdomain deflation • 1 vector represents left part of domain • 1 vector represents right part of domain • The eigenvector corresponding to the smallest eigenvalue is in the span of these two vectors. • Eigenvalues of and are almost the same, but the smallest is cancelled now. Deflated CG method

28. Deflated CG method Deflation Techniques Results for the simple problem: • Less iterates are needed • Residuals go faster to zero when • using deflation GOAL: IMPROVE CONVERGENCE BEHAVIOUR OF THE IBRAHYM MODEL or REDUCE NUMBER OF ITERATIONS

29. Deflated CG method Results for IBRAHYM

30. Results for IBRAHYM Process • First: Subdomain deflation while storingmatrix Z and AZ • Problem: Memory limiting for large areas. • Optimization concerning memory: • Claypackages  Layer based deflation (nice structure of Z) • Re-using vectors. Deflated CG method

31. Results for IBRAHYM Process • Problem: No gain of wall-clock times. • Optimization concerning wall-clock times: • ‘Storing’ AZ using one vector. • Pointer-vector instead of IF-loops. Deflated CG method

32. Results for IBRAHYM Results (small area) • Claypackages  Layer based deflation • For 3 small areas (7x7 km): • 1.489.600 cells Deflated CG method

33. Results for IBRAHYM Deflated CG method

34. Results for IBRAHYM Results (small area) REMEMBER: MODFLOW uses outer and inner iteration loops stop criteria inner loop: less iterations needed also: less fluctuation in solution when varying maximal number of inner iterations Deflated CG method

35. Results for IBRAHYM Results (small area) Layer baseddeflation also‘solves’ for thefaults. Deflated CG method

36. Results for IBRAHYM Results (bigger area) bigger area is 18 x 18 km = 9.849.600 cells • less iterations needed • wall-clock time gained • error smaller Error = inflow – outflow also: less fluctuation in solution when varying maximal number of inner iterations Deflated CG method

37. Deflated CG method Conclusions & Recommendations

38. Deflated CG method Conclusions & Recommendations Conclusions: • In general: deflation works for modelling groundwater flow. • Less iterations are required. • Deflation preconditioner makes solution more robust. • Wall clock times can be gained, but depends strongly in code used. Recommendations: • Deflation in horizontal direction. • Code should be further optimised for Fortran: • Minimizing memory, • IF-statements, • Using smart mappings. • Cluster multiple layers in one subdomain.

39. Deflated CG method QUESTIONS?

40. Deflated CG method Storing AZ

41. Deflated CG method Storing AZ

42. Deflated CG method Conclusions & Recommendations • Conclusions: • In general: deflation works for modelling groundwater flow. • Less iterations are required. • Deflation preconditioner makes solution more robust. • Wall clock times can be gained, but depends strongly in code used. • Recommendations: • Deflation in horizontal direction. • Code should be further optimised for Fortran: • Minimizing memory, • IF-statements, • Using smart mappings. • Cluster layers in one subdomain.