## Section: New Results

### Numerical models and simulations applied to hydrogeology

#### Reactive transport

Participants : Jocelyne Erhel, Souhila Sabit.

This work is done in the context of the MOMAS GNR ( 8.1.1 ) and the contract with Andra ( 7.1.1 ).

Reactive transport models are complex nonlinear Partial Differential Algebraic Equations (PDAE), coupling the transport engine with the reaction operator. We consider here chemical reactions at equilibrium. We have pursued our work on a global approach, based on a method of lines and a DAE solver [16] , [12] , [25] . We started also to analyze SNIA and SIA approaches, to find out conditions on the time step.

#### Flow and transport in highly heterogeneous porous medium

Participants : Jocelyne Erhel, Nadir Soualem, Julia Charrier, Géraldine Pichot.

This work is done in collaboration with A. Beaudoin, from the University of Le Havre (who moved to the University of Poitiers in September 2010), J.-R. de Dreuzy, from the department of Geosciences at the University of Rennes 1 (who is on leave for three years at UPC, Barcelona, Spain) and G. Pichot (who was at the University of Le Havre until september 2010). It is done in the context of the Micas project ( 8.1.2 ).

We have pursued our work for simulating flow and solute transport in 2D domains, where the permeability field is highly heterogeneous and is a random field We have worked on a new method to generate the permeability field based on spectral simulation. This method allows to improve the criteria required for the simulation of the permeability field. For a log-normal exponentially correlated field, in the two dimensional case, we have defined several conditions on the domain size and the mesh size to ensure a satisfying permeability field generation. We are currently working on establishing those criteria for 3D simulations. A paper is in preparation.

Stochastic computations are performed by using a Monte-Carlo method. Concerning the numerical analysis of the Monte-Carlo method, it has been done in the case of an isotropic and constant dispersion tensor, under stronger assomptions, proving the convergence of the method and yielding an upper bound for the errors comitted by estimating the spreading and the macro-dispersion of the solute. In the case of the spreading we have

where N is the number of realizations of the permeability field, M the number of particles for each realization of the permeability field, t the mesh of the time discretisation for the computation of the particles trajectories and h the mesh of the spatial discretization for the computation of the velocity field. A paper is in preparation.

Our on-going research concerns also the numerical study of convergence when we vary the numerical parameters in (1 ). The main goal is to provide optimal parameters for studying the macro-dispersion. Large scale simulations are used as a reference for tuning the parameters. We run 100 Monte-Carlo simulations which take each 3.2 hours with 128 processors on the IBM Power 6 at IDRIS. We observe numerically that convergence of Monte-Carlo iterations is quite fast; we study how the ergodic properties of the random permeability field can explain this behavior. A paper is in preparation.

We have also extended the transport model to include a hydrodynamic dispersion effect [11] . We study the different numerical approaches for dealing with discontinuities in the dispersion coefficient.

#### Flow in 3D networks of fractures

Participants : Jocelyne Erhel, Baptiste Poirriez, Géraldine Pichot.

This work is done in collaboration with J.-R. de Dreuzy, from the department of Geosciences at the University of Rennes 1 (who is on leave for three years at UPC, Barcelona, Spain) and G. Pichot (who was at the University of Le Havre until september 2010). It is done in the context of the Micas project ( 8.1.2 ).

We simulate flow within fractures lying in an impervious rock matrix. Discrete Fracture Networks are complex 3D structures with 2D domains intersecting each other. A first challenge comes with the meshing of such networks, where the mesh must be of good quality and must not contain too many cells. In a previous work, we have designed a method to generate a conforming mesh of good quality. In fractured media, flow is highly channelled in a small number of fractures. These fractures need to be finely meshed, while others can be coarsened. In order to reduce the number of cells, each fracture is meshed independently, resulting in a non-conforming mesh at the intersections. We addressed this difficulty within the numerical method by developing a Mortar method. For networks where intersections do not cross nor overlap, this Mortar method is based on pairwise relations [15] , [37] . The next step was to generalize the previous method to any stochastic network. This was a challenging task, since fractures can be highly intricated. We developed a new approach based on a combination of pairwise Mortar relations with additional relations for the overlapping part. A paper has been submitted. The Mortar method has been implemented in the H2OLAB platform in the module called MORTAR (see section 5.4 ).

Once the network is meshed and a mixed element method is applied, flow computations consists in solving a large sparse linear system. To test different solvers, we have developed an interface, called D3_ FLOW_ SOLVER (see section 5.4 ). A second challenge is to take advantage of the matrix structure to use domain decomposition methods and efficient parallel solving. Indeed, the matrix is naturally partitioned into blocks such that a Schur complement can be defined. Then the reduced system can be solved by a Preconditioned Conjugate Gradient method. An efficient method is the so-called Neumann Neumann preconditioner. However, floating subdomains can arise and lead to singular blocs. We have developed two strategies to overcome this difficulty. First, the subdomain decomposition consists in defining connected groups of fractures, so that the kernel of each block is at most of dimension 1. This is achieved by using a graph partitioning implemented in the Schotch library Scotch ). Second, we defined an algebraic update of the blocks associated with floating subdomains. Moreover, this approach allows to mutualize the factorizations of the blocks in both the preconditioning matrix and the Schur matrix [38] . We have developed a Matlab prototype implementing this method. It has been validated with a set of various fracture networks. We currently develop a C++ library called SOLVER_ SCHUR, integrated in the platform H2OLab (see section 5.4 ). Parallel computations are done with a AIX platform (IBM Power 6) at IDRIS.

#### Uncertainty Quantification methods

Participants : Julia Charrier, Jocelyne Erhel, Mestapha Oumouni.

This work is done in collaboration with A. Debussche, from ENS-Cachan-Rennes and Ipso INRIA team. It is done in the context of the Micas project ( 8.1.2 ). It is also done in collaboration with Z. Mghazli, from the university of Kenitra, Morocco, in the context of the Co-Advise and Hydromed projects ( 8.2.1 , 8.3.1 ).

In our applications described above, we use stochastic models and rely on uncertainty quantification methods. We have pursued the work on elliptic partial differential equations with random coefficients. We focused on the case of a lognormal homogeneous permeability field. This work applies in the particular case of an exponential covariance, which is one of the most frequently used model. We have then to deal with a permeability field whose trajectories do not have the usual regularity and are neither uniformly bounded from above nor from below.

We are in particular interested in numerical methods based on the approximation of the random permeability field using a truncated Karhunen-Loëve expansion. In the research rapport [42] , after proving that the solution belongs to L^{p}(, H1_{0}(D)) for any finite p , we provide both strong and weak error result estimates for the error on the solution resulting from the truncature. Moreover we give bounds for the spectral collocation error and the finite elements error. This work has been presented in several workshops and seminars: [35] , [17] , [20] , [19] , [18] , and has been pursued with an improved finite element error estimate.
This work is now extended to the numerical analysis of the multi-grid Monte Carlo method, in collaboration with Robert Scheichl, Aretha Theckentrup and Ivan Graham, from the university of Bath, GB. A paper is in preparation.

We have defined a numerical method based on the truncated Karhunen-Loëve expansion of the inverse of the random permeability field. This technique allows to approximate the mean of the solution by a projection. It is very efficient for 1D problems. A paper has been submitted.

#### Inverse problems in hydrogeology

Participants : Sinda Khalfallah, Jocelyne Erhel.

This work is done in collaboration with J.-R. de Dreuzy, from Geosciences Rennes, and A. ben Abda, from LAMSIN, Tunisia. It is done in the context of the Hydromed and Co-Advise projects ( 8.3.1 , 8.2.1 ).

We study two types of inverse problems in hydrogeology. The direct transient model is governed by classical flow equations and relates transmissivity with hydraulic head. We assume a constant known porosity.

The first type of problem is a so-called data completion problem, whith missing data on some part of the boundary and overdefined data on the other part. We have investigated methods based on an energy norm [13] . We have also defined a method based on a fictituous domain decomposition technique [33] , [39] , [28] .

The second type of problem concerns the identification of the transmissivity in saturated aquifers. We have studied a methodology based on pilot points.