Motivation
=========

A correct design and control of the irrigation cycles is very important for the correct and efficient grow of crops. It is well known that a correct design of irrigation cycles increases its efficiency and reduces the amount of water required. An important objective in the agricultural industry is to design an optimal irrigation plan depending on crop requirements. The irrigation cycles have an effect on the quality and quantity of water. It becomes more relevant in areas of the world where water is scarce. Nowadays, several scientists are involved in different research projects that focus on finding alternative water sources for irrigation. One idea is to reuse human waste water [1]. This concept is based on the pollutant removal effect of the sun, the soil and the plants itself. However, it requires an important control of the polluted water evolution in the subsurface. 

Numerical modelling is an important tool to help engineers and scientists design structures, simulate experiments, etc. In the context of this study, it can be used to design irrigation structures. It serves also to evaluate the propagation of pollutants or nutrients available in the irrigation water through the subsurface.

We present a useful app designed in Comsol that simulates the effects of different irrigation cycles on the soil saturation. It also can simulate the evolution of the pollutants dissolved in the water through the subsurface attending to its degradation. The application solves the Richard’s equation (eq 1) to obtain the groundwater flow field and its saturation. 

$$
\rho(\frac{C_m}{\rho g}+S_eS)\frac{\partial p}{\partial t}+\nabla \cdot \rho (-\frac{k_s}{\mu}\mathbf{k_r}(\nabla p+\rho g \nabla \mathbf{D}))=Q_m
$$

The advection diffusion equation (eq 2) is used to simulate the transport of pollutants and nutrients. 

$$
\phi \frac{\partial c_i}{\partial t}+\mathbf{u}c_i=\nabla \cdot [(\mathbf{D_{D,i}}+\mathbf{D_e})\nabla c_i]+R_i
$$

The model is also prepared to evaluate the evolution of the dissolved elements involved in the organic matter degradation chain (figure 1). 


<figure className="image">
 <StrapiImage
 src="/uploads/com_irri_app_fig_1_98f8084b34.png" className="mx-auto" alt=""
 />
 <figcaption>Figure 1. Sequence of important redox processes at pH 7. Figure taken from: Appelo and Postma, 2014 [2].</figcaption>
</figure>


The chemical reactions considered in the model are the following:

$$
O_2+CH_2O \leftrightarrow CO_2+H_2O
$$
$$
MnO_2 + 0.5CH_2O + 1.5CO_2 + 0.5H_2O \leftrightarrow Mn^{2+}+2HCO_3^-
$$
$$
Fe(OH)_3 + 0.25CH_2O + 1.75CO_2 \leftrightarrow Fe^{2+}+2HCO_3^-+0.75H_2O
$$
$$
SO_4 + 2CH_2O \leftrightarrow H_2S + 2HCO_3^-
$$
$$

This sequence is progressive which means that the denitrification, e.g., does not begin until all the oxigen is depleted. The reactions are included in the model as second order kinetics (Eq 4) where the rate is defined by the user.


\frac{\partial O_2}{\partial t}=-k_{O_2}[O_2][CH_2O]
$$
$$
\frac{\partial NO_3}{\partial t}=-k_{NO_3}[NO_3][CH_2O]
$$
$$
\frac{\partial MnO_2}{\partial t}=-k_{MnO_2}[MnO_2][CH_2O]
$$
$$
\frac{\partial SO_4}{\partial t}=-k_{SO_4}[SO_4][CH_2O]
$$
$$
\frac{\partial CH_2O}{\partial t}=\frac{\partial O_2}{\partial t}+0.5\frac{\partial NO_3}{\partial t}+0.5\frac{\partial MnO_2}{\partial t}+0.25\frac{\partial Fe(OH)_3}{\partial t}+2\frac{\partial SO_4}{\partial t}
$$

The app is fully parametrized. Regarding the geometry, it is up to the user to select the extension of the modelled area (Figure 2), the extension and location of the irrigation areas, etc. The app also considers different materials for the subsurface. Once the geometry is generated, the app allows the user to assign the properties of the materials such as the permeability or the porosity and introduce the irrigation planning which is read from an external file. Finally, the initial concentration of the different species is assigned together with the composition of the irrigation water. 

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_2_4d49ba9199.png" alt=""/>
 <figcaption>Figure 2. Example of a geometry generated in the app with one irrigation zone and three different materials. Red lines and points show parameters such as thickness or position that are set by the user.</figcaption>
</figure>

The app is prepared to show several types of plots such as 2D slices that show the evolution of the saturation, the lead concentration, etc. It also shows the evolution of the groundwater level and water saturation at given observation points.

Example of application
==================

This example considers two permeable materials separated by a small aquitard. The hydraulic conductivity of the materials can be observed in Figure 2. Two irrigation zones are included (Figure 3). The model considers the irrigation of the small crop only. It considers a regional groundwater flow from south to north. Under these boundary conditions, the computed groundwater flow streams from the irrigation zone located south-west to the north north-east (Figure 4 and Figure 5). The groundwater level is higher in the upper aquifer than in the lower. There is a general flow from top to bottom through the aquitard. However, this flow between aquifers is small compared to the horizontal flow. This lower flow retards the movement of the concentration through the aquitard ( Figure 6). Figure 3 to 8 show different types of plots that can be used to visualize the results in the application. In addition, the app also contains different tables to export the results in a quantitative way.


<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_3_062efd8df6.png" alt=""/>
 <figcaption>Figure 3. Finite element grid formed by tetrahedra and prisms.</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_4_5ae687cd9e.png" alt=""/>
 <figcaption>Figure 4. Computed hydraulic head.</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_5_3b9723faa3.png" alt=""/>
 <figcaption>Figure 5. Cut plane showing the hydraulic head and the flow direction.</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_6_95c22195a8.png" alt=""/>
 <figcaption>Figure 6. Chloride concentration in a cut plane. Arrows show the flow direction.</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_7_68be4bf646.png" alt=""/>
 <figcaption>Figure 7. Cl mol/m3.</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/com_irri_app_fig_8_94dff2d315.png" alt=""/>
 <figcaption>Figure 8. Isosurface Cl 1.2 mol/m3.</figcaption>
</figure>

The possibility of generated custom apps with Comsol allows the user to generate general models such as the one presented here that can be applied for different crops and areas. 

<figure className="image">
 <StrapiVideo width="725" height="395" src="/uploads/comsol_app_irrigation_video_ddd59299b5.mp4" type="video/mp4" />
 <figcaption></figcaption>
</figure>



References
==========

[1] *Hussain I.; L. Raschid; M. A. Hanjra; F. Marikar; W. van der Hoek*. 2002. **Wastewater use in agriculture: Review of impacts and methodological issues in valuing impacts**. Working Paper 37. Colombo, Sri Lanka: International Water Management Institute.

[2] *Appelo, C.A.J and Postma.D,* 2014. **Geochemistry, groundwater and pollution**. 2nd edition.

large_com_irri_app_fig_0.png

medium_com_irri_app_fig_0.png

small_com_irri_app_fig_0.png

thumbnail_com_irri_app_fig_0.png

com_irri_app_fig_0.png

We present a useful app designed in Comsol that simulates the effects of different irrigation cycles

A computational app for a proper evaluation of the irrigation effect over the aquifers

Model description
===============
 
A simple 2D domain has been used to explore the methodology for implementing the water age calculation.
 

Figure 1 shows the geometry, parameters and boundary conditions of the example. It can be seen that groundwater will enter the domain by the upper left part (left red segment) and discharge through the upper right part (right red segment).

 

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_1_6dbbff6e87.png" alt=""/>
 <figcaption>Figure 1. Sketch of the model.</figcaption>
</figure>



 
Mode 1: Transient simulation
=====================
Equation solved
------------------
 

For this model, a simple transport equation of an instant pulse of tracer (entering the domain) must be solved. Then, the equation that must be solved to computed groundwater age was proposed by Goode D.J.,1996:

 
$$
A=\frac{\int_{0}^\tau c(t) t dt}{\int_{0}^\tau c(t) dt}
$$


Where *A* is the “Water age”, *t* is the time, *C* is the concentration of tracer and $\tau$ is the total duration of the tracer pulse.

 
Implementation in Comsol
-----------------------------------
The transport boundary condition is a prescribed concentration on the recharge zone (hydraulic head=1), and open boundary in the discharge zone (hydraulic head=0). No flux is prescribed in the rest of the boundaries.

The required instant pulse of tracer was implemented by using the Comsol capability to create a time function associated to the boundary condition. This step function could be used as time function “step(t)” in the solute transport physics. A prescribed concentration boundary condition, associated to the time function “step(t)” was used to generate the pulse, as shown in Figure 2 and 3.

 <figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_2_df8d64e13a.png" alt=""/>
 <figcaption>Figure 2. Example of step function in Comsol.</figcaption>
</figure>

The model is then solved in a transient mode until all tracer concentration exits the model domain. 
 
 <figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_3_09a6ce429f.png" alt=""/>
 <figcaption>Figure 3. Implementing the step function in the transpor boundary conditions.</figcaption>
</figure>


The time integration equation is implemented in the results section using the function timeint(from,to,variable) where from and to are t=0 and the time when all tracer exits the domain, respectively. In this case the expression needed is:

 

**timeint(t0,T,c*t)/timeint(t0,T,c)**

 

Where t0 and *T* is the initial and final time respectively for the integral and *c* is the solute concentration.

The result of use this expression is a surface showing the computed water age, as shown in Figure 4.

 <figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_4_cb683ddc94.png" alt=""/>
 <figcaption>Figure 4. Results of transient method. Time in days.</figcaption>
</figure>
 

Mode 2: Steady state
==================
Equation solved
----------------
 

The equation solved for steady state method is:

 $$
\frac{\partial A}{\partial t}= 1 - \nabla \cdot A \frac{\mathbf{q}}{\theta} + \nabla \mathbf{D} \nabla A + \frac{F}{\theta \rho}
$$


Where A is the water age (seconds), q is the Darcy flux (m/s), D is the dispersion tensor (m2/s), (Theta) is the porosity, ρ is the density (kg/m3) and F is the source term (equal to 0 in the conceptual model). In steady state the temporal term dA/dt is equal to 0.

Implementation in Comsol
-------------------------------
 

In this method the transport boundary conditions are open boundary for zones with prescribed hydraulic head and no flux for the rest of boundaries. Also a source term in all domain is needed to obtain the water age, the value of the source term is equal to one in the equation, if we use Darcy flux (q) and effective diffusion coefficient the source term must be the porosity, as shown in Figure 5.

 
 <figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_5_64d35d5069.png" alt=""/>
 <figcaption>Figure 5. Solute source condition with a value equal to porosity.</figcaption>
</figure>




 

The computed water age is shown in Figure 6. It can be seen that computed results with the steady-state methodology is almost identical to what was computed at the transient approach (Figure 4).

 


 <figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/groundwater_age_figure_6_d2ad378405.png" alt=""/>
 <figcaption>Figure 6.Results of steady-state method. Time in days.</figcaption>
</figure>


 

Summary and conclusions
=====================

Two different approaches have been implemented in Comsol Multiphysics for calculating the water age in an aquifer. The steady-state methodology is easier, more accurate and much faster than the transient approach.

Acknowledgements 
===============

My more sincere thanks to Elena Abarca by her guidance and good advices and to Paolo Trinchero by his help to solve the steady state approach.

 

References
==========

*Goode D J,* (1996). **Direct simulation of groundwater age. **Water resources research, vol 32.n No.12, pages 289-296

medium_groundwater_age_figure_6.png

small_groundwater_age_figure_6.png

thumbnail_groundwater_age_figure_6.png

groundwater_age_figure_6.png

This post is intended to illustrate how to calculate the water age in a volume (3D model) or surface (2D model) based in the methodology proposed by Goode D.J. (1996).

Computing groundwater age in Comsol Multiphysics: implementation of transient and steady-state methodologies

Background
==============

An interface coupling two standalone codes, the geochemical calculator PHREEQC (Parkhurst and Appelo, 2013) and the multiphysics software COMSOL (COMSOL, 2013), has been developed recently by Amphos 21. The interface, called iCP (Nardi et al., 2013), is able to simulate complex thermo-hydro-mechanical-chemical (THMC) models.

The reactive transport process is solved using the operator-splitting (OS) technique (Strang, 1968) with a sequential non-iterative approach (SNIA). In this method the conservative solute transport is computed in a first step with COMSOL (coupled with other physical processes if desired) while the geochemistry is solved in a second step with PHREEQC.


<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/molal_solute_transport_figure_1_b02a96a138.png" alt=""/>
 <figcaption>Figure 1. Sequential non-iterative approach scheme of iCP, with the Molal solute transport interface highlighted.</figcaption>
</figure>

Despite COMSOL offers its own conservative solute transport interface, a new interface, “Molal solute transport”, have been developed for the specific needs of iCP. The interface is consistent with the units and the reactive transport formulation used in iCP. 

 Customizing COMSOL
==================
The interface was implemented using the Physic Builder toolkit of COMSOL. Within the Physic Builder environment, users can design customized interfaces through an interactive desktop and without the need of coding. The tailored interfaces can be exported as a plugin that can be added into the plugins folder of COMSOL and loaded every time the program starts. The interface will be, then, available as any other COMSOL module.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/molal_solute_transport_figure_2_5324c9711b.png" alt=""/>
 <figcaption>Figure 2. The Molal solute transport interface development in the Physics Builder environment.</figcaption>
</figure>

In the finite element method the governing equations have to be defined in weak form. Three different ways are available to introduce the weak form in COMSOL. In this case the weak formulation was implemented to allow a variable number of unknowns.

The Molal solute transport interface is able to model transport by advection and diffusion-dispersion. The interface can be used for steady-state and time-dependent analysis and in all kind of geometries supported by COMSOL.

The interface
=============
All features in the interface were tailored to meet the requirements of iCP. For instance, the term “Molal” comes from molality (mol/kg water), the default concentration units in PHREEQC; and the diffusion coefficient is constant for all the dependent variables, because in iCP formulation the transported elements are components and not “real” chemical species. 

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/molal_solute_transport_figure_3_c4e9c9609b.png" alt=""/>
 <figcaption>Figure 3. Example of a model implementation with the Molal solute transport. </figcaption>
</figure>


Three default nodes are added to the Model Builder when using the interface: Conservative Transport (the governing equation), No Flow (the default boundary condition), and Initial Values. Many other boundary conditions can be added: Fixed Concentration, Outflow, Open Boundary and Mass rate.

 The interface has been tested in various models, both reactive and conservative transport, performing similarly to other solute transport codes.

Acknowledgements
================

This work was made possible by the funding from the European Union's through the CO2-ReactMarie Curie Initial Training Program and has been conducted at Amphos 21 in close collaboration with Albert Nardi

References:
================
 

[1] *Parkhurst, D.L. & Appelo, C.A.J.,* 2013. **Description of input and examples for PHREEQC version 3—A computer program for speciation, batch-reaction, one-dimensional transport, and inverse geochemical calculations. ** U.S. Geological Survey Techniques and Methods, book 6, chap. A43, p.497 p. Available at: http://wwwbrr.cr.usgs.gov/projects/GWC_coupled/phreeqc/

 

[2] *COMSOL, 2013*. Comsol Multiphysics. Version 4.3a, Available at: www.comsol.com.

 

[3] *Nardi, A., Idiart, A., Trinchero, P., de Vries, L. M. and Molinero, J.* (2013). ** ‘Interface Comsol-PHREEQC (iCP), an efficient numerical framework for the solution of coupled multiphysics and geochemistry’, **Computers & Geosciences, submitted. Available at: http://image-modelling.net/icp.

 

[4] *Strang, G., *1968. **On the construction and comparison of difference schemes. ** SIAM J Numer Anal, 5, pp.506-17.

medium_molal_solute_transport_figure_3.png

small_molal_solute_transport_figure_3.png

thumbnail_molal_solute_transport_figure_3.png

molal_solute_transport_figure_3.png

In this post, Álvaro Sáinz, Marie Curie PhD student at Amphos 21, blogs about a COMSOL interface he has built to simulate conservative transport in porous media.

Molal solute transport: customized physics interface for iCP

The problem
================
Disposal in a deep underground repository is currently considered by several countries as the preferred option for long-term management of spent nuclear fuel. For regulatory reasons the safety of such a disposal facility must be demonstrated quantitatively through Performance Assessment calculations. One of the factors that have the potential to significantly influence the results of Performance Assessment is the radionuclide release rate from the waste immediately upon contact with groundwater (the so called **Instant Release Fraction**). Performance Assessment calculations indicate that instant release is important for the overall safety assessment of the facility, but there are significant uncertainties about the nature of this phenomenon and how it should be represented in Performance Assessment.

 With the above in mind, the **FIRST-Nuclides** Collaborative Project was funded by the European Commission under the 7th Framework Programme of the European Atomic Energy Community (EURATOM). The project has been initiated to provide better (conceptual and quantitative) understanding of the processes involved in instant release of radionuclides from spent nuclear fuel.

In the context of instant release, the knowledge of the rate of spent fuel saturation with water is key as the release of some radionuclides is contingent upon wetting of the fuel surfaces by water. Therefore, the instant release rate of these radionuclides depends on the rate of spent fuel saturation. An improved understanding of the rate of spent fuel saturation with water will help elucidate radionuclide instant release both under repository conditions and with regard to laboratory experiments conducted under near-room conditions.

Approach
=======

Under repository conditions the exposure of spent fuel to groundwater would be conditional upon the failure of multiple interconnected engineered barriers. The processes involved in such a phenomenon are multiple and coupled, making the system highly complex.

As a point of departure therefore, the saturation of a single idealised spent fuel pellet was considered (henceforth referred to as the *Reference Pellet*). This situation can be thought of as representative of batch laboratory leaching experiments performed on complete or fragmented fuel pellets, which are performed in order to measure radionuclide release rates from spent fuel upon contact with water. The results of modelling of saturation with water of a single spent fuel pellet are directly applicable to the interpretation of laboratory leaching results, but are also the necessary first step towards more complex multiphysics models capable of representing the fate of spent fuel upon wetting under repository conditions.


The conceptual model assumed that the pellet was composed of an impervious to water UO2 matrix and a network of (thermal) cracks, which can be filled with water (Figure 1). By saturation with water of the pellet is therefore understood the process of filling up of these cracks with water.

 

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/modelling_saturation_figure_1_9d08583700.png" alt=""/>
 <figcaption>Figure 1: The conceptualised Reference Pellet (left – lines indicate intersection of cracks with the pellet’s surface) and the finite element mesh (right) of the idealised crack system. Note that water flow is considered to occur through the cracks only (the fuel matrix is impervious to water). Dimensions shown are in metres.</figcaption>
</figure>

Implementation in Comsol
======================
Water penetration into the interconnected system of cracks was simulated by solving the Richard’s equation (in the pressure-based form). This was implemented in Comsol utilising the Fracture Flow interface (calculation of flow in discrete fractures), which was modified by manually implementing a storage-retention term and a set of constitutive equations (the water retention model) so as to account for variably-saturated flow within the cracks. More details on the conceptual and numerical details of the model can be found in [1].

The primary variable solved for is water pressure, which, through the water retention model, links directly with water saturation. In order to represent the physics of water penetration a negative water pressure (positive capillary suction pressure) was prescribed homogenously within the cracks, while a positive pressure (here corresponding to 5 cm of water column) was applied to all external crack edges. The driving force for water penetration into the cracks is water pressure (capillary suction) gradients (Figure 2), which disappear after the saturation process is complete.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/modelling_saturation_figure_2_070a9e6cf2.png" alt=""/>
 <figcaption>Figure 2: Successive stages of saturation with water (from left to right) of an interconnected system of cracks in a single spent fuel pellet. The colour scale represents effective saturation (-), which varies from the initial value of 0.1 ‰ (deep blue) to unity (red) at full saturation.</figcaption>
</figure>


Insights Gained
=============
The model enables quantification of the time needed to completely saturate a given system of cracks, which may vary depending on the statistical properties of the cracks (e.g. their number, geometry, interconnectivity and aperture), and depending on the initial and boundary conditions used (e.g. representing the presence of partial cladding or depth at which water penetration occurs). This is useful information as it makes it possible to compare the crack system surface wetting rate with the rate of radionuclide release observed in laboratory leaching experiments. Preliminary results of this on-going work indicate that, indeed, the time scales for cracks saturation with water are comparable with those operationally defined for instant/fast release of some radionuclides (hours to weeks).

We are currently working on an implementation of tracer release from the crack surface and transport by diffusion through the crack system towards an external reservoir (experimental batch solution or rock formation water in case of repository conditions). Also, the saturation model is being extended to include two distinct, but interconnected, populations of cracks (of different mean apertures). Ultimately, we intend to include relevant chemical reactions (e.g. dissolution of the fuel matrix at later stages) and upscale the model to the rod assembly scale. Work progress and results obtained in this work will continue to be presented and discussed at this blog.

Acknowledgements
================
The research leading to these results has received funding from the European Union's European Atomic Energy Community's (Euratom) Seventh Framework Programme FP7/2007-2011 under grant agreement n° 295722 (FIRST-Nuclides project). 

This work was carried out with important contributions from Amphos 21 colleagues: Andrés Idiart, Lara Duro and Olga Riba.

References
================
[1] *Pękala et al.* 2012. **Modelling of spent fuel saturation with water – approach, preliminajry results and potential implications. **In: KIT Scientific Reports 7639, 1st Annual Workshop Proceedings of the Collaborative Project “Fast/Instant Release of Safety Relevant Radionuclides from Spent Nuclear Fuel (7th EC FP CP FIRST-Nuclides).

medium_modelling_saturation_figure_2.png

small_modelling_saturation_figure_2.png

thumbnail_modelling_saturation_figure_2.png

modelling_saturation_figure_2.png

In this post, Marek Pękala describes a recent application of Comsol at Amphos 21 to study the saturation with water of spent nuclear fuel.

Modelling saturation with water of spent nuclear fuel in the context of its geological disposal

The basics
=============

Helium originates mainly from three sources: the atmosphere, the mantle and as a result of radioactive decay processes taking place in the crust. Each source of helium is characterized by a particular 3He/4He ratio. The isotopic composition of helium in the atmosphere is characterized by a 3He/4He ratio of about 1.4x10-6. Radiogenic helium is primarily 4He, with a 3He/4He ratio in the terrestrial rocks closer to 10-8, and primordial helium in the deep earth is enriched in 3He with a variable isotopic composition (from 1.1x10-5 in the upper mantle to 8.7x10-5 in the lower one). Helium in the crust originates from the alpha decay of natural uranium and thorium series, whereas 3He is formed by thermal neutron capture by 6Li.

Helium production in the crust is dominated by the alpha decay of the uranium and thorium decay chains, and therefore directly proportional to the concentration of 235,238U and 232Th in the rocks.

The Helium production rate in the rock can be evaluated as:

$$
G=\frac{\rho}{V_m}\left( 1.19 \cdot 10 ^{-16} [U] +2.88 \cdot 10 ^{-17} [Th] \right)
$$

where r is the rock bulk density (kg/m3), Vm the molar volume of the gas (m3/mol) and [U] and [Th] are the uranium and thorium contents of the rock matrix (ppm), respectively. Then, if we know (measure) the concentration of dissolved helium in a groundwater sample, we can estimate how much time this water had (theoretically) been in contact with the rock to reach such Helium content. But this is assuming all the Helium comes from radiogenic origin without mantelic and/or atmospheric contribution. And also assuming a closed 0-dimensional system, i.e. no diffusive transport of Helium occurs in the rock porewater.

New analytical solution for interpretation of He data in deep groundwater
============================================================

In our paper we present the details of a new analytical solution for modeling Helium production and diffusion in a 1D finite profile. It is important to remark the finite nature of the model domain because this is one of the key innovations. So far the only available analytical solution, known as the “Andrews model” (Andrews 1985), was for an infinite thickness of the crust. The other key innovation with respect to “Andrew’s model” is the possibility of discretizing the crust in different domains with heterogeneous [U] and [Th] compositions (i.e. to include different rocks in the model). The analytical solution looks like this:

$$
c\left( z, t_r \right) = \frac{2d}{D_p} \sum^{\infty}_{k=1} \frac{ sin \left(\alpha_k \frac{z}{d} \right)}{\alpha_k^2} \left[ 1 - exp \left( \alpha_k^2 \frac{D_p}{d^2} t_r \right) \right]
\times \\
\left\{ \frac{G_1 d}{\alpha_k \phi} \left[ cos \left( \alpha_k \frac{d_l}{d} \right) - cos \left( \alpha_k \frac{d_u}{d} \right) \right] + \left( -1 \right)^k j_0 \right\}
$$

Where $D_p$ is the pore diffusion coefficient, tr is the residence time of groundwater, $G$ is the Helium production, f is the rock porosity, $J_0$ is the mass flux of Helium coming from the mantle (bottom boundary condition) and d is the total depth of the model (i.e. the deep boundary of the finite domain where $J_0$is applied as a boundary condition). In the example we have published, the discretization of the domain is done by 3 intervals like:

$$
 G \left( z \right) =
 \begin{cases}
 0, & \ 0>z>-d_u \\
 G_ 1, & -d_u \geq z>-d_l \\
 0, & -d_l \geq z >-d
 \end{cases}
$$

where $d_u$ is the depth of the “shallow” aquifer with rapid groundwater renewal (and we assume the production of he is negligible), and dl is the maximum depth of the rock unit containing [U] and [Th]. It is worth noting that we adopt this “discretization” to mimic the application example shown below, but it can be changed or extended to any other settings.

 

Applying the new methodology at Olkiluoto (Finland)
=======================================

Olkiluoto is the site where POSIVA (the finish company for nuclear spent fuel management) is conducting field investigations to evaluate the suitability of constructing an underground facility to dispose high-level radioactive waste. Olkiluoto is located on the coast of Gulf of Bothnia in south-western Finland. The two main lithological units of the Olkiluoto site are (1) high grade metamorphic rocks and (2) igneous rocks. Hydrogeological and geochemical investigations shows that the site has a “shallow” part (from surface to about 250 meters) with active hydrogeological circulation, and a deeper zone with groundwater movement is much lower and limited to some water conducting fracture zones.

Radiogenic helium production at Olkiluoto is related with the presence of granites in the bedrock. Recent seismic studies have suggested that in the Olkiluoto area, different granites extend to a depth of at least 4 km. Thus, we have considered two models, one where a constant production rate is applied from the upper boundary (which is located at 250 m depth, after the zone with more active hydrogeological circulation) down to 4000 m depth, and the second where in situ production is present only in the first 1000 m of bedrock. A deep source of helium ($J_0$), which mimics deep-seated granite production, is located at the bottom of the domain.

The Uranium and Thorium contents of the rock are well known, as well as diffusion coefficients and porosities. Then, we did a fitting exercise by calibrating the groundwater residence time, tr and the deep source of helium, $J_0$. Both models were able to fit very well the Helium concentrations measured in the field, as shown in the next figure (dots are field measurements and lines are the calibrated models).

 
<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/00_Groundwater_Age_Helium_2_0694b5b652.png" alt=""/>
 <figcaption>Figure 1. Comparison among measured Helium concentrations and the best fit calibration with the new analytical solution. The fitting parameters are the residence time of groundwater and the Helium flux at the bottom boundary.</figcaption>
</figure>

The resulting fitted parameters are the same for 2 models, which indicates that extending the Helium production deeper than 1 kilometer is not too sensitive. The best fit was obtained for a tr (residence time) of 3.2x107y, and $J_0$(Helium flux from the mantle) of 3.1x10-9 mol/m2y.

So far we do not have data about Helium isotopic ratios in the groundwater. Future determinations of such isotopes could help us to quantify the mantelic contribution and then test the validity (or at least the degree of confidence) about our estimation of 30 million years of residence time for the groundwater in the granitic rocks of Olkiluoto.

ACKNOWLEDGEMENTS
=====================

I would like to acknowledge to all my coauthors of this paper. It was a great project and I think we all enjoy very much working together. Paolo Trinchero was the leader of the project at that time, together with Anne Delos who worked very hard in the first stages. Marco Dentz was the one really behind the maths for the analytical solution and Petteri Pitkanen provided the best expertise for Olkiluoto geological settings and the behavior of dissolved gasses in the groundwater. Ignasi Puigdomenech (SKB) was also a key person for this project, and provided high level comments and constructive suggestions to the first versions of the report and paper.

References
=========

*Trinchero, P.; Delos, A.; Molinero, J.; Dentz, M. & Pitkanen, P. *(2014). **Understanding and modelling dissolved gas transport in the bedrockof three Fennoscandian sites.** Journal of Hydrology, 512, 506–517.

*Andrews, J. 1*(985). **The isotopic composition of radiogenic helium and its use to study groundwater movement in confined aquifers.** Chemical Geology 49, 339–351.

medium_00_Groundwater_Age_Helium_1.png

small_00_Groundwater_Age_Helium_1.png

thumbnail_00_Groundwater_Age_Helium_1.png

00_Groundwater_Age_Helium_1.png

Dating groundwater is something that hydrogeologists like to do for different practical purposes

How old is groundwater? Estimation of residence time with Helium

Background
=========
Groundwater inflows into tunnels may constitute a potential hazard as well as an important factor controlling the rate of advancement in the excavation. In fact, one of the major practical difficulty during tunnel construction is often associated with the presence of groundwater. Predictions of groundwater inflows are also needed to design suitable drainage systems. Furthermore, drawdowns produced by tunnel construction may induce land subsidence, water table decline and environmental impacts on rivers and wetlands.

In the literature, there are several analytical expressions to calculate water discharges into tunnels but a major development was presented recently by Perrochet and Dematteis (2007). The program TunFlow makes use of the above-mentioned analytical expression to compute groundwater infiltration into tunnels and galleries. The program allows an unlimited number of geological zones (heterogeneous system) to be accounted for. Furthermore, fully transient hydrogeological conditions (accounting for the velocity of tunnel advancement) are considered. TunFlow is easy to use and is especially suited for the evaluation of groundwater filtrations in scoping studies for tunnels construction.

 

Input parameters
=============
The analytical expression developed by Perrochet and DeMatteis (2007) to calculate the groundwater tunnel inflows and used in TunFlow is shown below:

$$
Q(t) = 2 \pi \sum^{N}_{i=1} H(t-t_i) \int^{v_i\left( t-t_i \right)}_{0}
\frac{K_i s_i H \left( L_i - x \right)}
{ln \left( 1 + \sqrt{ \frac{\pi K_i}{S_i r_i^2} \left( t-t_i - \frac{x}{v_i} \right) } \right) } dx
$$

The parameters of this formula are:

- $t$ [d]: time
- $L$ [m]: length of the section of a given geological layer
- $v$ [m/d]: drilling velocity
- $K$ [m/d]: permeability of the considered layer
- $S$ [1/m] Specific storage of the considered layer
- $r$ [m]: tunnel radius
- $s$ [m] level of water table (distance from the tunnel roof)
 

Before running the model, it is preferable to create an Excel sheet with all the information requested by the code. The numbers of layers is specified in the Graphical User Interface (GUI) of TunFlow (1 and 2, **Figure 1**).The number of time steps must be also selected (3, **Figure1**) . By default, 1000 time steps are used. It is worthwhile noting that increasing the number of time steps improves the accuracy of the solution but also increases the computational costs of the simulations (i.e. simulation time).

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/tunflow_1_03_72cf2245f5.png" alt=""/>
 <figcaption>Figure 1. Screenshot of the Graphical User Interface (GUI) of TunFlow.</figcaption>
</figure>

Once all the parameters are imported in TunFlow, the program can be executed (4, **Figure 1**).

 

Results and output options 
=============
 

Tunflow provides diferent outputs. Calculated data can be exported in txt and xlsx file formats (1 and 2, **Figure 2**).On the other hand, the program also plots the results as the tunnel inflows Q (l/s) versus time (d) (3, **Figure 2**).

 

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/tunflow_1_04_d4f5cb40c1.png" alt=""/>
 <figcaption>Figure 2. Results and output options of TunFlow.
</figcaption>
</figure>


The input parameters can be also exported as an xlsx file (1 and 2 in **Figure 3**) and used for new simulation runs.

 


<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/tunflow_1_05_55b81820f6.png" alt=""/>
 <figcaption>Figure 3. Preparation of the input parameters in an Excel spreadsheet.
</figcaption>
</figure>

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/tunflow_1_06_20968f3d42.png" alt=""/>
 <figcaption>Figure 4. Import of the previously created Excel file into Tun Flow.
</figcaption>
</figure>

 
Acknowledgements
===============
This work has been conducted at Amphos 21 in collaboration with Santos Jiménez, Paolo Trinchero and Jorge Molinero.

References
=========
*Perrochet, P., Dematteis, A.*, 2007. **Modeling Transient Discharge into a Tunnel Drilled in a Heterogeneous Formation.** Ground Water 45 (6), 786–790.

medium_tunflow_1_01.png

small_tunflow_1_01.png

thumbnail_tunflow_1_01.png

tunflow_1_01.png

In this post, Jordi Font, Senior Consultant and Project Manager at Amphos 21, blogs about a program (developed at AMPHOS 21)

Tunflow 1.0 - How to calculate groundwater inflow into tunnels and galleries

Background and motivation
======================

In performance assessment (PA) studies of nuclear waste repositories, numerical simulations are used to study the potential impact of the release and transport of radionuclides in the near and far field. These simulations are typically carried out using a lumped approach: the different and highly non-linear underlying radionuclide retention processes are represented by a single parameter, the distribution coefficient (Kd), which relates the radionuclide mass retained in the solid phase to its aqueous concentration.

Kd-based simulations, which are particularly appealing due to their numerical efficiency, are mathematically robust but usually rely on an important assumption: the Kd value of each radionuclide is constant in time. This assumption could be critical when geochemical conditions are changing dramatically (e.g. in the transition from temperate to glacial conditions), as distribution coefficients depend on the pH, redox conditions and major chemistry of the system.

We have recently developed and tested a framework (interface to Marfa - iFM) that constitutes a significant improvement over the traditional constant Kd-based models. The basic premise of this approach is that the “background” geochemistry (excluding radionuclides) is first be pre-calculated using a reactive transport code in order to generate a dataset of Kd values in space and time that are representative of the actual geochemical conditions in the considered spatial and temporal domain. This information is then used by the computer code Migration Analysis of Radionuclides in the Far Field (MARFA) [1,2] to calculate radionuclide transport based on these “intelligent” Kd values.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/smart_Kd_figure_1_bb6a5d7027.png" alt=""/>
 <figcaption>Sketch showing a typical radionuclide transport problem. Radionuclides, which are released at repository depth, move along a number of transport pathways and are discharged into the biosphere.</figcaption>
</figure>

The tool: iFM version 2.0
=====================
The figure below shows a workflow of iFM v2.0.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/smart_Kd_figure_2_2eeb82baaa.png" alt=""/>
 <figcaption>Functionality overview of iFM v2.0. </figcaption>
</figure>

The numerical tool consists of two sets of numerical calculations:

- FASTREACT [3] simulates transport and geochemical reactions along the same set of transport pathwasy used by MARFA. This simulation only accounts for the background geochemistry (without radionuclides).
- PHREEQC batch calculations. At each time step, the background geochemical composition computed along the pathways (previous step) and an arbitrarily low concentration of the considered radionuclide(s) are used to perform zero-dimensional calculations. These batch calculations use a mechanistic description of the relevant sorption processes (e.g. surface complexation, cation exchange, etc) and provide as output quantitative information about the distribution coefficients, Kd, of each radionuclide at each location and at the considered simulation times.
 

Testing the tool – a synthetic example
======================

To test the robustness of the tool, we have solved a “typical" radionuclide waste performance assessment (PA) problem: the simulation of the migration of radionuclides from a deep geological repository to the surface.

A snapshot of a groundwater velocity field has been taken from an existing regional-scale three-dimensional hydrogeological model. Using particle tracking simulations, a number of transport pathways have been delineated.

Radionuclides have been released during a small time interval compared to the time frame of the problem. Initially, glacial conditions are assumed to prevail in the whole domain. After around 20,000 years, a different and less diluted water, representative of new temperate conditions, enters the model domain from the release boundary and propagates undisturbed along the streamlines.

The numerical exercises have been solved using (1) iFM+MARFA (see Figure below) and (2) FASTREACT.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/smart_Kd_figure_3_11bb3789da.png" alt=""/>
 <figcaption>Sketch that shows the approach and the different processes included in the iFM+MARFA simulations.</figcaption>
</figure>



The strontium breakthrough curves computed using iFM+MARFA and FASTREACT are shown in the Figure below. When strontium is overtaken by the geochemical perturbation, it is mobilized and moves nearly unretarded. Both simulations capture this effect and their results closely match each other. Thus, the verification exercise can be considered successful.

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/smart_Kd_figure_4_671d85dc4e.png" alt=""/>
 <figcaption>Calculated breakthrough curve for strontium (iFM+MARFA: blue solid line; FASTREACT: orange dashed line).</figcaption>
</figure>

Acknowledgements
===============
This work is a collaborative project between SKB and POSIVA, the organizations in charge of managing nuclear fuel and radioactive waste in Sweden and Finland, respectively. Amphos 21 and Scott Painter want to thank Lasse Koskinen (POSIVA) and Jan-Olof Selroos (SKB) for their support and very helpful comments.

Paolo Trinchero also wants to thank the coauthors of this work: Luis Manuel de Vries, Jorge Molinero and Scott Painter.

References
=========
[1] *Painter, S., V. Cvetkovic, J.-O. Selroos, and O. Pensado.* (2007). **Time domain particle tracking methods for simulating transport with retention and first-order transformation. **Water Resources Research 43(9), W01406.

[2] *Painter, S., and Mancillas, J. *(2013). **MARFA User’s manual: migration analysis of radionuclides in the far field, **Working Report 2013-01, Posiva Oy.

[3]* Trinchero, P., Molinero, J., Román-Ross, G.* (2013) **A streamline-based approach for the solution of multicomponent reactive transport problems,** SKB report R-10-45.

medium_smart_Kd_figure_1.png

small_smart_Kd_figure_1.png

thumbnail_smart_Kd_figure_1.png

smart_Kd_figure_1.png

Distribution coefficients, also known as Kd values, are lumped parameters that contain information about the radionuclude-soil/rock interaction processes

An “intelligent” Kd approach for the simulation of radionuclide transport in the geosphere - iFM version 2.0

<div className="">
 <a className="inline-flex items-center text-gray-500 bg-transparent border border-solid border-gray-500 hover:bg-gray-100 hover:text-white active:bg-gray-600 font-bold uppercase text-sm px-6 py-3 rounded-full outline-none focus:outline-none mr-1 mb-1" href="https://nitrate.amphos21.com/" target="_blank">Launch Experiment</a>
</div>



Motivation
=============
Understanding the behaviour of pollutants in aquifers is a common issue we have to face as hydrogeologists. Numerical models are useful tools that allow us to predict the displacement and evolution of the pollutants in the subsurface. For this modelling task we usually use well-known software, but we find important **keeping simple approaches in mind since they can be useful for learning about numerics insights and physics understanding**. 

To this end, we have developed a simple python program to compute the Advection-Diffusion equation (Bear, 1972) with first order decay and degradation chain [Eq.1] in 1 dimension using finite differences. This **has been applied to simulate the nitrification chain**, the degradation of ammonium generates nitrite and this specie when degrades produces nitrate. ($NH_4^+ \rightarrow NO_2^- \rightarrow NO_3^-$).

Description
=============
The system modelled is represented by the following partial differential equations:

$$
\tag*{[Eq. 1]}
\begin{cases}
\frac{\partial NH_4^+}{\partial t}=D\frac{\partial^2 NH_4^+}{\partial x^2}-\frac{\partial}{\partial x }(aNH_4^+)-\lambda_1 NH_4^+\\
\frac{\partial NO_2^-}{\partial t}=D\frac{\partial^2 NO_2^-}{\partial x^2}-\frac{\partial}{\partial x }(aNO_2^-)+\lambda_1 NH_4^+-\lambda_2 N0_2^-\\
\frac{\partial NO_3^-}{\partial t}=D\frac{\partial^2 NO_3^-}{\partial x^2}-\frac{\partial}{\partial x }(aNO_3^-)+\lambda_2 N0_2^--\lambda_3 N0_3^-\\
\end{cases}
$$

Here, $a$ is the fluid velocity (water in the studied case) (m/s), $D$ is the diffusion coefficient ($m^2/s$) and $\lambda_i$ is the half-time of the $i$ species where 1 refers to $NH_4^+$, 2 is $NO_2^-$ and 3 is $NO_3^+$.


Obtained results (Figure 1) have been validated by comparing with the analytical solutions reported by (Cho et al. 1971).

<figure className="image">
 <StrapiImage className="mx-auto" src="/uploads/numerical_vs_analytical_282d15042c.png" alt=""/>
 <figcaption>Figure 1. Analytical solution from (Cho et al. 1971) (solid lines) vs numerical solution (dashed lines).</figcaption>
</figure>


Numerical solution
==============
Previous to the **finite differences** scheme development, it is useful to define some terms:

- $L$ = Length of the domain.
- $N$ = number of grid points.
- $· x = L/N$ spatial discretization.
- $i$ is the position at the domain being, therefore, $i*· x$ where $i \in [O,N-1]$ is the real position.
- $T$ = Total simulation time.
- $k$ refers to a determined timestep being $· t$= time step and $k*·t$ the real time.

As can be seen in [Eq. 1] the proposed equation has a diffusive and a convective term. Using a Crank-Nicolson scheme for the diffusive term:

$$
	\frac{-D}{2· x^2}u^{k+1}_{i-1}+(\frac{1}{· t}+\frac{D}{· x^2})u_i^{k+1}-\frac{D}{2 · x^2}u^{k+1}_{i+1}=\frac{D}{2· x^2}u^{k}_{i-1}+(\frac{1}{· t}-\frac{D}{· x^2})u_i^{k}+\frac{D}{2 · x^2}u^{k}_{i+1}
$$

Here, $k$ is the current time and $i$ is the cell in the domain.

And for the advection term:
$$
	\frac{\partial au}{\partial x} \approx \frac{a}{4· x}(U_{i+1}^k-U_{i-1}^k+U_{i+1}^{k+1}-U_{i+1}^{k+1})
$$

Combining this two expressions, we obtain the advection diffussion equation in a cell point (i) at a time (k):

$$
 	\frac{U_i^{k+1}-U_i^k}{· t}=\frac{D}{2· x^2}(U_{i+1}^k-2U_i^k+U_{i-1}^k+U_{i+1}^{k+1}-2U_i^{k+1}+U_{i-1}^{k+1})+ \\ 
\frac{a}{4· x}(U_{i+1}^k-U_{i-1}^k+U_{i+1}^{k+1}-U_{i-1}^{k+1})-\lambda U_i^k
$$

We define $\sigma = \frac{D· t}{2· x^2}$ and $\rho = \frac{a· t}{4 · x}$. Reordering the previous approximation:

$$
 (-\sigma +\rho)U_{i-1}^{k+1}+(1+2\sigma)_i^{k+1}-(\sigma +\rho)_{i+1}^{k+1}=(\sigma -\rho)U_{i-1}^k+(1-2\sigma)U_i^k+(\sigma + \rho)U^k_{i+1}+· t (-\lambda U^k_i),\\
 i=1,..,N-2.	
$$

This equation holds for all grid points containing in the domain $i\in [1,N-2]$. The expression in the borders (0 and N-1) is given by the boundary conditions. In this case, we assume Neumann conditions in both boundaries.

The resulting expression for $i$ equal to $0$ and $N-1$ is:

$$
	 j=0: (1+\sigma +\rho)U_0^{k+1}-(\sigma+\rho)U_i^{k+1}=(1-\sigma -\rho)U_0^{k}+(\sigma+\rho)U_i^{k+1}+· t \lambda(U_0^{k})
$$

$$
 	 j=N-1: (-\sigma +\rho)U_{N-2}^{k+1}+(1+\sigma -\rho)U_{N-1}^{k+1}=(+\sigma -\rho)U_{N-2}^{k}+(1-\sigma +\rho)U_{N-1}^{k}+· t \lambda(U_{N-1}^{k})
$$

 The previous expression can be written as a linear system in matrix notation:

$$
 	AU^{k+1}=BU^{k}-f^k
$$

 Where:
$$

 A[NxN]=
 \begin{bmatrix}
 1+\sigma +\rho & -\sigma+\rho & 0 & 0 & 0 & 0 & 0 & 0\\
 -\sigma-\rho & 1+2\sigma & -\sigma+\rho & 0 & 0 & 0 & 0 & 0\\
0 & -\sigma-\rho & 1+2\sigma & -\sigma+\rho & 0 & 0 & 0 & 0\\
 ... & ... & ... & ... & ... & ... & ... &...\\
 ... & ... & ... & ... & ... & ... & ... & ...\\
 0 & 0 & 0 & 0 & -\sigma-\rho & 1+2\sigma & -\sigma+\rho &0\\
 0 & 0 & 0& 0 & 0 & -\sigma-\rho& 1+2\sigma& -\sigma+\rho \\
 0 & 0 & 0 & 0 & 0 & 0 & -\sigma-\rho & 1+\sigma-\rho
 \end{bmatrix}
$$

$$
 B[NxN]=
 \begin{bmatrix}
 1-\sigma -\rho & \sigma-\rho & 0 & 0 & 0 & 0 & 0 & 0\\
 \sigma+\rho & 1-2\sigma & \sigma-\rho & 0 & 0 & 0 & 0 & 0\\
 0 & \sigma+\rho& 1-2\sigma & \sigma-\rho & 0 & 0 & 0 & 0\\
... & ... & ... & ...& ... & ... & ... & ...\\
 ...& ... & ...& ... &... & ... & ... & ...\\
 0 & 0 & 0 & 0 & \sigma+\rho & 1-2\sigma & \sigma-\rho &0\\
 0 & 0 & 0& 0 & 0 & \sigma+\rho & 1-2\sigma& \sigma-\rho\\
 0 & 0 & 0 & 0 & 0 & 0 & \sigma+\rho & 1-\sigma+\rho
 \end{bmatrix}
$$

 $U^{k+1}$ and $U^k$ are the solution at a given time k and k+1.

$f$ is the source term which in our case is a vector of length $N$ with the following expression:

$$
f[N]= U^k\begin{bmatrix}
 	-· t \lambda ,-· t \lambda,....,-· t \lambda
 \end{bmatrix}
$$


Conclusions
==============
The proposed numerical scheme and its implementation are able to reproduce the results of the analytical solution reported by (Cho et al. 1971). The model describes the dissolved species movement through advection, diffusion and takes into account the degradation by first order decay.

References
===============
- *Cho C. M.*, 1971. **Convective transport of ammonium with nitrification in soil.Canadian journal of soil science.**
- *Bear J.,* 1972. **Fluid dynamics in porous media.**

medium_numerical_vs_analytical.png

small_numerical_vs_analytical.png

thumbnail_numerical_vs_analytical.png

numerical_vs_analytical.png

The current experiment presents a Crank-Nicolson scheme to solve the advection-diffusion equations

Modelling Advection-Dispersion with Decay Chains

[Amphos 21](https://www.amphos21.com) is thrilled to announce an exclusive 1.5-hour course on [GibbsStudio](https://gibbsstudio.io/), the geochemical modeling framework that's transforming how we integrate Phreeqc into our workflows.


<figure className="image">
 <StrapiImage
 src="/uploads/gibbsstudio_dark_mode_dea5dacd6a.png" className="mx-auto" alt="GibbsStudio Dark mode"
 />
</figure>

 
📅 **Mark Your Calendars!**
- **Date**: 29th of May
- **Duration**: 1.5 Hours
- **Pricing**: Free
 
🔍 **What to Expect?**

This engaging course is designed to provide you with a comprehensive understanding of GibbsStudio's capabilities. This session is specially tailored for those who regularly utilize Phreeqc and wish to enhance their workflow with GibbsStudio's advanced features.
 
🚀 **Course Highlights:**

- Hands-on experience in GibbsStudio's environment.
- In-depth guide through key modules: Phreeqc, Akva, Predominance, among others.
- Real project examples developed at Amphos 21, showcasing applications in various real-world projects.
- Live Q&A session to address your questions and curiosities.
 
🎯** Who Should Attend?**

Geochemists and hydrologists who regularly work with Phreeqc.
 
💡 **Why Join?**

Seize the opportunity to boost your Phreeqc productivity with GibbsStudio, acquiring practical knowledge and skills that will empower your research or professional endeavors in geochemical modelling.
 
📌 **Reserve Your Spot!**

Spaces are limited. Ensure you secure your participation by fill in [the form](/form/gibbsstudio-course-registration-request/).
 
Stay tuned for more details, and feel free to reach out to us at techlabs@amphos21.com for any queries or further information.

large_gs_logo.png

medium_gs_logo.png

small_gs_logo.png

thumbnail_gs_logo.png

gs_logo.png

Barcelona Science and Amphos 21 are thrilled to announce an exclusive 1.5-hour course on GibbsStudio.

Hydrogeology