Deflated Conjugate Gradient Method for modelling Groundwater Flow

Deflated Conjugate Gradient Method for modelling Groundwater Flow in a Layered Grid Lennart Ros Deltares & TU Delft datum 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)

Outline Deflated CG method 12/12/2021 2

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 Deflated CG method 12/12/2021 3

Introduction Deflated CG method 12/12/2021 4

Introduction Deltares January 1 st 2008 Deflated CG method 12/12/2021 5

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

Introduction Geohydrology § The driving force for groundwater flow is the difference in height and pressure. § To represent this difference we introduce the concept of hydraulic heads, h [L]. Deflated CG method 12/12/2021 7

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

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

Introduction MODFLOW: § MODFLOW is a software package which calculates hydraulic heads. § Developed by the U. S. Geological Survey. § Open-source code: everyone can use and improve this program § Rectangular grid and uses cellcentered variables. § Quasi-3 D model. Deflated CG method 12/12/2021 10

Introduction IBRAHYM: § groundwater model developed for several waterboards in Limburg. § 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 of solver. Deflated CG method 12/12/2021 11

Equation, Discretization & Method Deflated CG method 12/12/2021 12

Equation, Discretization & Method Governing Equation: Where: Deflated CG method hydraulic conductivities along x, y, and z coordinate axes [LT-1], h head [L], W volumetric flux per unit volume representing sources and sinks of water [T-1], Ss specific storage of porous material [L-1], t Time [T] 12/12/2021 13

Equation, Discretization & Method Finite Volume Discretization: Deflated CG method 12/12/2021 14

Equation, Discretization & Method Discretized Equation Using Finite Volume Method: Where: Deflated CG method 12/12/2021 15

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

A Simple Example Deflated CG method 12/12/2021 17

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

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

A Simple Example Observations for a simple testcase in Matlab: Deflated CG method 12/12/2021 20

A Simple Example Observations for a simple testcase in Matlab: Smallest eigenvalue: 0. 00010283296716 Next eigenvalue: 0. 04870854847951 Deflated CG method 12/12/2021 21

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 Deflated CG method 12/12/2021 22

Deflation Techniques Deflated CG method 12/12/2021 23

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

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 Deflated CG method 12/12/2021 25

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 Deflated CG method Subdomains 12/12/2021 26

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 12/12/2021 27

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 Deflated CG method 12/12/2021 28

Results for IBRAHYM Deflated CG method 12/12/2021 29

Results for IBRAHYM Process § First: Subdomain deflation while storing matrix 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 12/12/2021 30

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 12/12/2021 31

Results for IBRAHYM Results (small area) § Claypackages Layer based deflation § For 3 small areas (7 x 7 km): § 1. 489. 600 cells Area original code deflation code 1 (yellow) 62 53 2 (blue) 164 121 3 (orange) 287 168 Deflated CG method 12/12/2021 32

Results for IBRAHYM Deflated CG method 12/12/2021 33

Results for IBRAHYM Results (small area) REMEMBER: MODFLOW uses outer and inner iteration loops Concerning maximal number of inner iterations: stop criteria inner loop: less iterations needed also: less fluctuation in solution when varying maximal number of inner iterations Deflated CG method 12/12/2021 34 Max. inner iteration Original code Deflation code 20 566 271 30 383 199 40 335 179 50 287 168 75 230 153 total number of iterations when varying the maximal number of inner iterations per inner loop

Results for IBRAHYM Results (small area) Concerning faults: Faults Original code Deflation code With 287 171 Without 254 168 Layer based deflation also ‘solves’ for the faults. Total number of iterations when we take and do not take faults into account Deflated CG method 12/12/2021 35

Results for IBRAHYM Results (bigger area) bigger area is 18 x 18 km = 9. 849. 600 cells Original code Deflation code § less iterations needed § wall-clock time gained § error smaller Iterations 1082 757 Wall-clock time 1664. 40 sec 1542. 20 sec Error -0. 31188 E-02 0. 49375 E-03 Error = inflow – outflow also: less fluctuation in solution when varying maximal number of inner iterations Deflated CG method 12/12/2021 36

Conclusions & Recommendations Deflated CG method 12/12/2021 37

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. Deflated CG method 12/12/2021 38

QUESTIONS? Deflated CG method 12/12/2021 39

Storing AZ Deflated CG method 12/12/2021 40

Storing AZ Deflated CG method 12/12/2021 41

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. Deflated CG method 12/12/2021 42